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RESPONSE OF AN ELASTIC DISK TO IMPACT AND 
MOVING LOADS*+ 


By A. CEMAL ERINGEN? (Purdue University, Division of 


Engineering Sciences, Lafayette, Indiana) 
Received 20 December 1954] 


SUMMARY 
With the use of Fourier transforms a class of elasto-dynamic problems concerning 
disks have been solved. The disk is subjected to various types of dynamic loadings 
at the rim. The cases of impact and moving loads are studied in detail. 


1, Introduction 

CYLINDRICAL roller bearings in high speed mechanisms are subject to 
dynamic loadings. Yet the usual design procedure is based on Hertz 
formulae which are the result of elasto-static considerations. In many 
other instances, gears, rollers, or disks are subject to impact or moving 
loads. If we neglect the coriolis terms we can also bring the rolling disks 
in contact into the category of disks subject to moving loads. Thus the 
aim of the present paper is to obtain the solution to this class of elasto- 
dynamic problems concerning the disk. The dynamic load is applied to 
the rim of the disk. Two normal concentrated dynamic loads at the two 
ends of a diameter moving or otherwise are special cases. 

Some solutions of free oscillation of cylinders have been known since 
the time of Pochhammer (1), and later Pickett (2), J. Mindlin (3), T. Ghosh 
4). Similarly, the problem of rotating disks has attracted attention of 
many authors (see, for instance, Lamb and Southwell (5), Timoshenko 
and Goodier (6), Love (7)). It seems, however, that the forced oscillation 
problems concerning disks and cylinders have escaped attention, excepting 
apaper by J. Mindlin (8), which consists of generalities in a problem related 
to ours. 


The present method is applicable to ring problems and to the plane 


elasto-dynamic problems concerning circular holes. These problems will 
be treated in later papers. 


t The research contained in this paper was obtained during investigations for a research 
project sponsored by the Office of Naval Research. 


- > . . ‘ 
} Associate Professor of Division of Engineering Sciences. 


(Quart. Journ. Mech. and Applied Math., Vol. VIII, Pt. 4 (1955)] 
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2. Formulation of the problem 
The equations of motion of plane homogeneous isotropic media in terms 
of plane polar coordinates 7 and @ and the time ¢ are (7) (see Fig. 1) 
yuy (A+ 2pu)A .—2ur-w 9 
Wy = (A+2pu)r—1A 9+2pw, 


where u(r, 6,¢) and v(r,6,t) are components of the displacement vector 


and « are the Lamé constants, and y is the mass density per unit volum 





Fic. 1. Cireular disk. 


Dilatation A and rotation w are related to u, v by 
rA (ru) +9, 2rw (rv) .— Ug. 


Subscripts after a comma represent differentiation, i.e. cu cb, et 
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Elimination of uw and v between (1) and (2) leads to 
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These are the equations of dilatational and rotational waves. 


Components o,,., o,9, gg of the stress tensor are given by: 


7 AA +-2ut », om pr lug tr er(v/r)» ) ; 
t 
on AA+ 2ur-!(v 9g+u) 


The problem is to solve (1), (2) under a given o,, and o,9 at the rim 
of the disk: 
o,,(a, @, t) 0, (9, t), a,9(a, 4, t) T,) (0, t) 


subject to the condition that these surface tractions are in equilibrium 
each instant. 
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} The solution 
The peri dic solution of (3) with respect to 6 is 
\ S 1,, sinné+A, cos n6)Z,(p,) | 
(6) 
SN (B,. cos nd B,, sin n®)Z.y(p) | 
bs ( D ( ; D ) (p;), Pj x; TI (j 1, 2), 
re Z he cylinder function (9): A,,, B,,, C;,, and D,, are constants 
tegratic nd the barred quantities represent the Fourier transforms, 
| | e' F(t) dt. (7) 
( rmula to s 
, | P _ 
e~ "7 F(z) dr. (3) 
Subs 6) into (2) after taking the Fourier transforms of (2) and 
g the resulting part differential equations we obtain 
> * sin n@ l% (r)cos né | 
(9) 
> | cos no I3,,(7)sin né | 
l beg ' Z (py) B,,, 2nps*Z (ps) | 
, (10) 
*Z,(p1)+ Bj, 2p2 * Z)(p2) 
e p epresent differentiation. 
( 1a) j 6 ind (9 we obtain 
Z » { \ b, Ny (rz) |sin n6 
A,,, N,,(’7)+ Bs N3,,(77) |cos né 
~ S { » 3 b N,. (v7) leos n@ 
« tails hina (11) 
Ay, S;,(77)+ Bg,, So,,(77) sin nO 
Z ») b,,, Ts rz) |sin né 
A,,, T,,,(r7)+ By, Ty, (7) |cos né 
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where 
Min(77) = (A/24)Z,,(p1)+-(1—n? py *)Z,,(p1) + pr * Zp (p4) 
Np, (17) = 2npz? Z,,(p2)—2npz ' Z;,( po) 
S,,(77) = npy* Z,(p;)—npx * Z,(p) 
Son (77) = (1—2n* py *)Z,,(p2)+2p2 ' Z)(p2) is 
T,,,(77) (A/2u)Z,,(p1) +07 * Z,,(p1)—pr |Z; (1) 
(r7) N,,, (77) 





4. The dynamic tractions applied to the rim of a disk 


In the case of a disk, the stress and deformation components must ly 


finite at r = 0. Hence D,, = 0. Without loss of generality we also tak 
Ci, 1. This means that in all our formulae we must replace Z, by J 
We use boundary conditions (5) to determine the constants A;, and B 


Fourier’s theorem thus leads to 


Ain [ 2D, (ar )][ 83, (AT )Go. N3,, (7) Foe 








B,., | 2D, (a7) | S,,, (47) Gos , Nj, (a7) Foe 
Ag, = [ 24D, (a7) J" San (47) 600+ Non (47 Fog . 
(ie 
B... [: 2uD, ( ar )] 1 S, (AT) Ge Nin (A7)F os | (n ty ae 
A jo 3] A, j = 4 0 Bio 4 B;, |, 0 
D, (ar) = N,,,(a7)S2,,(a7) —Ng,, (at) S,,,(a7) 
where as 
Gy. = 7! | &(8,7)cos nO dé, 
0 
Gog = 7! | 6,(80,7)sin nO dé, 
0 
and 7,, and 7), are similarly defined. 
Various special cases are of interest: 
(a) Zero surface shear: 7). = 7, = 9. (14 
(b) Normal traction with central symmetry and (a): &). = 9. (15 
(c) (b) with constant amplitude over 0 <0 <a, te. 
‘ . 0<8@ x 
(15) and 6,(@,7) G(T) 
|\aw—a <0 7 
0 a<0< 7—x«. (Il 


We get Foe = 2Go(7)(7n) | 1+ (—)"Jsin na. (Ii 
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1) Impact load and (c): 


With (16) and lim 26,(7)aa P,(7) 


e get Do 2P,(7)/7a. (18) 


Impulsive concentrated load and (d): 


With (18) and P(t) P.d(t), where 8(t) is the Dirac delta function 





fined by . 
t 0 iz 
ee oe 8(t) dt = 1, 
10 t Q, , 
e find P(r) i Coe 2P,/7a, (19) 


| where P, is the amplitude of the concentrated load. 


cn 


Moving loads 
) Moving normal and lange ntial stresses: 


Moving loads can be represented by 





o, (6, t) o,(6— Qt), T,)(0, t) 7)(8—Q, t). (20) 
Fourier transforms of these are 

6,7) Jet tig (7 / Cd), 7) (0,7) Zeit Mg (7/Q,), (21) 

l ore ; 
(7/Q | ao,(dje-7?2 dd (in general) 
18) 
ic 
7 a,(d)cos(th/Q) db when o,(—¢) o9(¢), (22) 
ere s,(7 Q,) follows from s,(7/Q) by replacing 2 and o,(¢) by 2, and 


T)\O) respec tively. Hence 


(2 |rQ2i)[ (r/Q)?—n2]-1(e?#7#— 1) 54(7/Q) (23) 
. (2. 
In a Q)| (4 Q)?—n?| L(g27TiQ 1)s (7 Q)) | 


| y \< 
| 


and 7,. are obtained from (23) by writing Q, and s,(7/Q,) 


| The quantities + 


I place of © and s,(7/ 02) respectively. 


(b Me ving periodic loads 
v0 


(dt) (Q,./ 7a > p, cos n(@— Qt) +-q,, sin n(@— Qt) | 


n 


>, 7, cos n(@—Q, t)+-s sin n(6—Q, t) | 
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we substitute in the ty 


After taking Fourier transforms of o, and 7». 


equations following (13) to obtain &, This gives 


tu 


AG oe @ (p,,—tq, )d(QQqn — 7) +-(p,, +2, )d(—22n —7) 
AG» , (ip, J, )0( Sn 7) ( ip, q,,)d( On i . 
The quantities 7), and 7, follow from (25) by writing 7,, and s,, in place 


p, and q,. In obtaining (25) we have used the formal relation? (10) 
270(u) |< tu dt. 


now be obtained } 


Components of displacement and stress tensor cai 
combining (9), (11), (13), and (25) and taking inverse Fourier transforms 
Thus 





(270 v9 yu S Th (7rQn | p,, COS n(@— Qt) qd, sin n(0 Qr)| 
n—0 
u\*)(rQn)!| ie sin n(@ Q, f)--s,, cos n(A—Q, t)] 
(27apu,Qyrje = > cP (rQn)| —p, sin n(8— QU -+-¢q,, cos n(8— QE) | 
} 0 
Uy (rQn)| — nO £2, ¢) s sin n(@ 2, t)| 
we C0 )Fry > oy}, (rQn)| PP, COS n(@— Qt) Tn sin n(@ {2t) | 
7 0 ‘ 
of) (7Q, n)| —r,, sin n(@ Q, t)-+s,, cos n(A—Q, t)] 
(7a Yo)o,6 z TY(rQn)| p, sinn(@— Qt) +-q,, cosn(@— Qt) 
, 0 
7\?)(rQn)| r,, cos n(O—Q, t) s, sin n(O Q,1)] 
(74, Qo)o9a = X orn (rQn)| p,, cos n(8— QF) +g, sin n(A— Q4)] 
n 0 
of) (rQn)| r sin n(@ (), t) s, cos n(O QO.) 
where 
uD (rQn) [| D, (aQn) |] (ay rQn) 1S) (og rOn)S,, (aQn) 


2n(agrQn) SF, (ay rQn)S,, (aQn) | 
WMV Qn) = [| D,(aQn)]Y n(ay rQn) SF, (ay rQn)S,,, (aQn) 


? 


2(a,7Qn) 1S) (ay rQn)S,,(aQn)| ) (2 





a) (rQn) | D (aQn )| . N,,(rQn Sp, (a Qn) Ny, (rQy )S,,(aQn )] 
TY (rQn) [ D, (aQn)] ALS), (rQn)S,, (aQn)— Sy, (rQn)S,,, (aQn)| 
oy) (rQn) D (aQn)| WT, (rQn)Sy,(aQn) T,, (rQH)S,, (aQn) | 


This is justified in the sense of a distribution function (11). 











nd Ww 


nence 


whel 


aT) 


Hen 
men 
The 





| 











>} 


SS componel 











ELASTIC DISK TO [IMPACT AND MOVING LOADS © 391 


(Qn A Qn)S, (aQn)—N, (aQn)S, (aQn). (28) 


ind o};) are obtained from the corresponding 


</ 





script (1) by replacing Q by Q, and S,,(aQn) and S,,(aQn) 
und NV, (aQ respectively, except in D (aQn), where we 
cally Silt moving, CONCE nitrated loads (Fig. 2), 
which is one of technical importance, we have 7, — 0, 
0. The ncentrated loads can formally be represented by 
QO. za Ot) +-8(7—04- Qt) | 
QV, 7a) > p, cos n(@—Qt)+-q, sinn(@—Qt), (29) 
an ee 
4 
Q, 4 
2. Moving radial load. 
umplitude of each of the concentrated radial loads. From 
coefficients p,, and q,, are 
(2/7 tor nm even, “ 
p qd, Q. (30) 
: iQ) for n odd 
ether with (29) and (30), andr, 8, 0, gives the displace 


ts. Below we give displacement components. 


ined in an obvious manner. 
T=ApL/ Vy S Us, (7O2n)cos n(G@— Qt) 
Pape Q, SY oP(rQn)sin n(@—Qt) | 


| moving shear loads and other types of load combination 


obtained from (24) and (26). 
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6. Computation and discussion 
Computations have been carried out to determine the roots of 
D (aQn) = 0 


given by (28). The roots of this equation give the resonance speed for the 
concentrated load. 


The ratio of the resonance speed C, = aQ of the concentrated load to 


the dilatational wave velocity C, = a,;1 has been solved from D, = 0 for 
a steel cylinder with EZ = 30 10° p.s.i. and v = 0-3. The first four roots 
of D, = 0 for n 1, 2, 3, 4 are listed in the table given below.+ 


The computation was carried out on I.B.M. Card Program Calculator 
and are correct up to two decimal places. 


The values of nC,/C, = nax,Q (obtained by calculating the roots of 
D (aQn) = 0). Steel: E = 30x 108, vp = 0-3 


n D, = 0 D,=0 D,=0 D,= 0 
I I°5I 1°26 1°94 2°52 
2 3°43 2°35 3°21 4°05 
3 3°30 4°22 4°95 5°07 
4 5°33 5°06 6°22 7°19 


From this table it is seen that the smallest C,/C, is obtained for n = 2. 
This value is C,/C, = 1-26/2 = 0-63. We therefore expect a critical speed 
for the moving load in the neighbourhood of 0-63C, which will create 
resonance in the disk. 

The Rayleigh surface wave velocity for v = 0-29 is 0-9258 times the 
velocity of shear waves or 0-503 times the velocity of dilational waves 
Thus 0-63C, represents a velocity between the Rayleigh surface wave 
velocity and the shear wave velocity. Further we notice that this velocity 
is minimum for n = 2 rather than n 1. Examining (31) we find that 
we have only even terms. Hence diametrically opposite moving loads 
give a smaller critical speed than a single load. In the former case the 
fundamental mode n = 0 represents a uniform lateral extension which is 
not dependent on time. Thus it is static in nature and has no resonance 
frequency associated with it. The second mode is of cos 20 type and gives 
the minimum critical speed mentioned above. It is interesting to note 
that the minimum critical speed in the case of a flat semi-infinite plate (the 
Rayleigh surface wave velocity) is less than that of the curved surface 
(the disk in the present case). With this point in mind perhaps one cal 
classify the curved surfaces from the elasto-dynamic point of view, thoug! 
this process will, no doubt, present mathematical difficulties. 


+ For n = 0 we take the limit of the functions in (27) and find finite values. 
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The present results, of course, must be accepted tentatively until 


verified experimentally. 
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THE EFFECT OF VISCOSITY UPON THE CRITICAL 
FLOW OF A LIQUID THROUGH A CONSTRICTION 


By A. M. BINNIE (Triniiy College. Cambridge) 
Received 24 February 1955] 


SUMMARY 


The motions of a viscous liquid under gravity over a broad-crested weir, and wit} 
swirl under pressure through a convergent-divergent nozzle, are analysed on thi 
supposition that the whole of the discharge passes through a laminar boundary layer, 
The constrictions are assumed to be gradual; and the pressure gradients in the 
direction of streaming being favourable, the velocity distributions over the cross 
sections are taken to be of polynomial forms. When the momentum integra 
equations are satisfied, first-order differential equations are obtained which shoy 
that, as with inviscid flow, the motion is ‘critical’ in the sense that arbitrary inlet 
conditions may not be specified. The flow must be determined from the conditions 
that prevail at the critical cross-sections, which are displaced downstream from thi 
throats of the constrictions. The equations are difficult to solve because at tl 
critical cross-sections, where their step-by-step integration has to be begun, they 
take the indeterminate form 0/0, but a numerical example of each type of motior 
is given. It is proved that a long wave of small amplitude can maintain itself 
stationary on the moving liquid at the critical cross-sections. A comparison is mad 


with inviscid flow through both kinds of constriction. 


1. Introduction 

AMonG the well-known examples of the ‘critical’ flow of a liquid are flow 
under gravity over a broad-crested weir and swirling flow under pressure 
through a convergent-divergent nozzle. These motions have been examined 
by several investigators, and a theoretical review was given by Binnie (1 
who showed that, in general, the flow sets itself in such a way that the 
streaming velocity is equal at the throat to the velocity of a long wave 
This analysis is confined to an inviscid liquid; and, as in Reynolds’s classical 
treatment of gas flow through a nozzle, the constriction is supposed to be 
gradual so that transverse components of velocity and acceleration maj 
be neglected. 

The present paper is an inquiry into the consequences of assuming that 
the liquid is viscous instead of perfect and that the whole of the flow passes 
through a laminar boundary layer, the inlet length being ignored. The 
working is approximate in that, following Pohlhausen, we take the velocity 
distributions in the boundary layers to be of reasonable polynomial forms, 
and employ von Karman’s method of satisfying the momentum integrals 

The two cases mentioned above are examined in detail. Although it 1s 
of no practical importance, the broad-crested weir is worth considering 


(Quart. Journ. Mech. and Applied Math., Vol. VIII, Pt. 4 (1955)] 
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because the analysis is relatively simple; thus the basic ideas may readily 
be understood and numerical results obtained with comparative ease. For 
swirling flow through a nozzle, Taylor (2) made it clear that in many prac 
tical cases the boundary layer is so thick that most of the discharge passes 
through it. This conclusion was supported by Binnie and Harris (3), who 
nalysed the motion for an inviscid liquid and found that, if the swirl were 
sufficiently great, the discharge was reduced to zero, all the energy in the 
supply being used up in producing tangential velocity. It is a fact of 
omm servation that the discharge does not cease in this way, hence 
nder these conditions boundary layer flow becomes predominant, as in 
some of the experiments on gravity flow through a trumpet described by 
Binnie and Hookings (4). The method of solving the weir problem is there 
fore applied to the nozzle, but the working is much more complicated, for 
two momentum integrals must be satisfied, yielding a pair of simultaneous 
lifferential equations. In both cases an outline of the theory for frictionless 

on is included, partly for purposes of comparison and partly to demon 


strate that at the throats the derivatives of the velocities and the thick 


esses of the streams take the indeterminate form 0/0. 


2. Frictionless flow over a broad-crested weir 
We consider unit width of the broad-crested weir outlined in side elevation 


F 


g. 1. The liquid passes steadily from a large reservoir, the surface of 











H 
= 
I |. Inviscid flow over weir. 
which stands at a level H above a horizontal datum, over a weir (shown 
ded), and falls away freely on the downstream side. At a typical cross 
section, the depth AB of the stream is denoted by 6 and the height of 
the weir above the datum by hk. Horizontal distances are given by the 
coordinate The vertical component of velocity being ignored, the hori 
ontal ve ty wis easily shown to be uniform over A B, and its magnitude 


ined from Bernoulli's equation 


(2.1) 
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in which g is the acceleration due to gravity. The equation of continuity 
” Q = ud, (2.2) 
where Q is the discharge per unit width. On differentiating (2.1) and 
(2.2) with respect to x, we obtain 

udu dd. dh 

gdx | dx! dx’ 


nseed 0 dd 5 du ( ) 4) 
é “Vapi ; 2. 
dx dx 
When du/dz is eliminated from (2.3) and (2.4), the result is 
dh 
dd dx . 
— (2 )) 
dx ;_ 
gd 
53 dh 
dx i 
D3" 2.0) 
Ps ” " 
q 
We are concerned here with the type of motion in which the liquid is acceler- 


ated throughout its passage over the weir, therefore d3/dx is negative 
everywhere. In the approach to the throat dh/dx is positive, and in that 
region u?/(g5) must be less than unity. At the throat defined by dh/dx = 0, 
ds/dx ~ 0, therefore from (2.5) 

u? = gy, (2.7) 
in which the suffix ¢ denotes throat values. It follows from (2.1) and 


9 7 ‘ 
(2.7) that H- h, 38). (2.8) 


Then, if 6, and the variation of h with x are specified, the profile of the 
stream can be determined from the relation 


h,+-38,—h eS, (2.9) 


which is obtained from (2.1), (2.2), (2.7), and (2.8). 
Another, but more cumbersome, way of evaluating 5 is also possible. 

If 5, = (Q?/g)! is specified, we could integrate (2.6) commencing at the 

throat, but there d5/dx takes the form 0/0. On applying to (2.6) the usual 

method of dealing with this condition, we find that for this section 
5S S 2 1 

dé 7 8d h > (2.10) 

dx 3 da? 


and the integration of (2.6) may be begun with this value. 
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If from (2.3) and (2.4) we remove dé/dx, we obtain 


dh 
du u dx 
——-, (2.11) 
dx 6 us 
l 
go 


ind this expression also takes the form 0/0 at the throat. 


3. Viscous flow over a broad-crested weir 


The entire thickness of the stream is now taken to be influenced by 


viscosity. We employ the same notation as before, adding (as indicated 


Y 
dx | 
ez 









Fic. 2. Viscous flow over weir. 


in Fig. 2) wu, to denote the horizontal velocity at the free surface and y to 
denote vertical distance above the constriction. It will be assumed that 
the distribution of velocity over A B has the same form wherever AB may 
be situated; and, writing » for y6, we may suppose wu to be given by 

u u,(27 27° } ”*), (3.4) 


for this relation satisfies the conditions 


u 0 at n 0) 
u u, at I ‘ 
1 / (3.2) 
ai 
0 at ”) l 
y 
Then 
Q u dy xu,5, where a = 0:7; (3.3) 
| w*>dy = Buys, where B = 0-5825; (3.4) 
0 


where y 
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Alternatively, in place of (3.1), we may adopt the simpler form 
u u,(2H ”*), (3.6 

which has the same shape as the full line in Fig. 4. This relation conforms 
to (3.2) and yields 

. 0-6. B 0-53. y =. (3.7 
There is thus no great difference between the two sets of numerical 
coefficients. 

The horizontal momentum integral can be obtained from first principles 
by considering the horizontal force exerted on the element A BCD shown 
in Fig. 2, The density, viscosity, and kinematic viscosity of the liquid 
being p, p, and v, this force, reckoned positive in the direction of streaming, 
has three components: 

(i) the difference between the hydrostatic forces on A B and DC, which 
is —gpd (db dx) dx; 
(ii) the horizontal component of the normal pressure exerted by BC, 
which is —gpd(dh/dx) dx: 
(iii) the viscous force on BC, which is — pl eu cy|, » da. 


Therefore the momentum integral is 


fe) 
a dd dh ou 
u? dy go go v ‘ (3.8 
C2 dx da cy 
0 ES 
which with the aid of (3.4) and (3.5) becomes 
soi (u28) dd dh uy a4 
6 uzo —Jo go yy « (3.! 
Cx , dx dx ra) 
From (3.9) uw, is removed by means of (3.3), and we obtain 
dé S3(dh/dx)+ T aT 
; (9.10 
dx o°— S 
i 
: bQ- ne vyY 
where y= e 1 tt I (3.11 
xg XY 


This equation for determining 6 should be compared with (2.6), which is 
the corresponding result for inviscid flow. 

Again dé/dx is negative over the whole of the weir. Now since dh dx is 
positive on the upstream side of the geometrical throat and negative on 
the downstream side, there is a certain cross-section where the numerato! 
of (3.10) changes sign. At this section, which may be termed the effective 
throat, the denominator must change sign also; or, in other words, dé da 


takes the form 0/0. Enforcing this condition we have 


x5 dh 


AY 


—T, and & S. (3.12 
dx 
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. second result can be transformed to 


go 
1 >? 


p 


(3.13) 


ich is similar to (2.7), and again a definite relation exists between the 


i th and the surface velocity at the throat. From (3.12) it appears that 
dh T ies 
vos, nis (3.14) 
da S BQ (Bqo*) 


Thus the position of the effective throat, which is determined by (3.14), is 
ifted from the geometrical throat in the opposite direction to that of the 
sternal foree due to viscosity. This result agrees with the analysis, given 

Binnie (1), of the swirling flow of an inviscid liquid passing under 


ivity down a vertical trumpet; in that case the external force lay in the 


lirection streaming. and the effective throat was found to be moved 
ipwards, although the displacement was very small in the numerical 
examples considered by Binnie and Hookings (4). 

The integration of (3.10) is begun at the effective throat, and again, 
efore a ut can be made, the difficulty due to the 00 form must be 


ercome. On applying the usual method to (3.10) we obtain 


» 2 1a\- dh dd Oo dh ‘ - 
| , ‘). (3.15) 

La dadx’ 3 dz? 
Since dh dx and d*h dx® are negative, this quadratic equation possesses 


one ne itive root and no ambiguity arises. 
If 8 is eliminated instead of u, from (3.9), we arrive at an expression 

luid vhich can be changed to 
/ uy 53 (dh dx) - 


/ e) >? NS 


(3.16) 


[his bears the same resemblance to (3.10) as (2.11) does to (2.5). and it 
takes the form 0/0 at the effective throat. 
With the same assumptions and methods as before, it may readily be 
own that (3.13) also gives the relation between the depth and surface 
locity of a horizontal stream on which a long wave of infinitesimal 
nplitude can maintain itself stationary. We take 6 (now constant) to be 
»depth of the undisturbed stream, and consider a vertical cross-section 
where w’ is the small additional surface velocity due to a wave of elevation e. 


\cross this section the velocity distribution is assumed to be given by 
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in accordance with (3.1), therefore from (3.3) the discharge is 


d+e 
Q= | udy = a(u,+u’)(6+). (3.18 
0 
Similarly, for the undisturbed stream 
4 ‘ y 2(v . (¥ a . 
Ui Z——Zi— — di U,o, 3. g 
qd | | 5 5 5 | y Xs € ( l 
0 
5 
hence uty = bo nl (3.20 
o-+e€ 


In applying (3.8) to this problem, we notice that dh/dx is now zero, and we 

neglect the last term, supposing that it is balanced by some external force 

which maintains the steady motion. Using (3.17), (3.4), and (3.20), we have 
d+e 


[ u* dy = B(u,+u')*(6+e) = Bui a(1 — ‘}: (3.21) 
0 7 
ie “ie . de 
thus (3.8) reduces to — But —gb —, (3.22 
dx dx 
9 go 
or Uz . 3.23 
1 B ( 


in agreement with (3.13). 


4. Numerical example of flow over a broad-crested weir 


The weir is taken to be of symmetrical parabolic shape defined by 


h 10-4(200a—<a?). (4.1 
Thus it extends from « = 0 tox = 200 with a maximum rise equal to unit) 
The critical depth is assumed to be 0-1 and to occur at x = 135; and there 
dh/dx, which is given by 
dh _ 3x10 1(100—2), (4.2 
dx 


has the value —0-007. It follows from (3.12) that 7’ must then be 7 x 10° 


If the foregoing figures are in inches, this requirement is met by water at 


about 17° C. for which v = 1-660 10-3 in.2/sec. For from (3.12) and 
(3.14) 7 can be written in the form 
53\4 
T=y sr (4.3 
"\Bg 


and on inserting into this the above values of 5 and v together with 
B = 0-5825 and y = 2 given in (3.4) and (3.5), we obtain the result aimed 
at. Accordingly, with the aid of (4.2), (3.10) becomes 
» —4 O— r)532§-- Tv 6 
2x 10-4(101 ad? + r x 10-6 (4.4 


dd 
-§8—0-001 


da 
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chich from the initial values 6 0-1 in., x 135 in., is to be integrated 
hoth directions as far as 2 0 and x 200 in. The solution has very 
andly been obtained for me by Dr. J. C. P. Miller, whose result is given in 


the satisfactory profile of the stream shown in column 5. It will be seen 
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umn 3 of Table I. The addition of the corresponding values of h yields 


it at small values of x, in the region most remote from the start of the 
lculation, the surface is very slowly falling; and, as it should, the slope 
ntinuously increases 

To afford a comparison, an example of frictionless flow has been added 
Table 1 by means of (2.9). In this the value of 6, at x 100 has been so 
osen that. at a 135 is the same as in the viscous case. From the 
sults set out in columns 4 and 6 it appears that, as expected, the profile 
s steeper, and the effect of viscosity in holding up the stream is well 
uked at the downstream end of the weir. For the first two entries in 
lumn 3 differ by 0-1901 and those in column 4 by 0-1907, whereas the 


st two entries differ by 0-00258 and 0-00444 respectively. 


TABLE I 
De pths and profile Ss (in inches) of VISCOUS 


and inviscid streams passing over weir 





ry — 
l i} 
I I 4 I 56 
1°O77 I I 7 
9005 I sO I »s 
a 2 ! c 2 > 
) I if I 
50910 1°203 I 
11 5 1*20074 is 
1°196 I 
7 I'l 55 I 
75 1°1754!1 1*20975 
027 I°r54! 1°15027 
7 14970 1°12157 I*13970 
1°0750 T*O5005 
7 10705 [°O1457 1°O1 795 
°) 10000 0°97750 0°97750 
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(3.1) may be obtained. The equation of motion corresponding to (3.8) js 


cu ae C2u 


u g _-(6+h)+v (4.5 
Cx dx 


ey? 
At the solid boundary (3.1) makes 6?u/éy? equal to zero, whereas the right 
hand side of (4.5) should vanish. To meet this boundary condition as we 
as those given by (3.2) the profile 

u u,(An- By? -C'n) ${ 
was tried at x 130. At this section u, 7-788 in./see. from (3.3) an 
(3.13), and with the aid of the third-difference formula given in /nii 


yielded « = 0-6725, B = 0-5404, y = 2-0705. These results should } 
compared with (3.3)—-(3.5). The effect of imposing the correct bounda 
condition is small at this section, but the final entries in column 5 sho 
large changes in 6-+-/, and there the discrepancy is greater. In that regior 
, 

of the kinetic energy to be concentrated near the free surface, but i 
inviscid flow the whole cross-section of the stream has equal acceleratior 
and with viscosity in operation the same tendency is to be expected. T 
redistribution by viscous action is not immediate as is assumed in (3.1). bi 
the error becomes smaller as the state of equilibrium between the gravit 
and the viscous forces is more nearly attained. Column 3 shows that 

the downstream end of the constriction the alterations in 6, and cons 
quently in uw, are proportionately less than elsewhere. The motion 

tending towards steady flow down an inclined plane, for which the exa 


velocity distribution is known to be given by (3.6). 


nozzle 
We consider steady flow from a reservoir maintained at a pressure head / 
which is sufficiently high for gravity effects to be negligible, through t! 
circular nozzle shown in Fig. 3. The tangential velocity v at any radius 
is related to the swirl constant 2 by the equation 
Q) 
- 
At a typical cross-section, the radii of the nozzle and of the core 4 
denoted by a and b, and the radial thickness 8 of the stream is given | 


From the table some estimate of the errors involved in the approximatio; 


polation and Allied Tables (6) the derivative of 5+-A was found to ly 
~0-006980. After A, B, and C had been evaluated, the calculations 


however, another source of error is least. Equation (3.1) supposes most 
| 


5. Frictionless swirling flow through a convergent-divergen! 
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2 Q 
» a is convenient to write v (5.3) 
N | h 
the tangential velocity at the core. To obtain non-dimensional results 
take Oo 
R 9 (5.4) 
9 - 
’ ()> 
kK ¢ : (5.5) 
- aa 
1 
D 
SS Uhl Viihihe 1117 T TTI Fry 
$ Oo / ie J 
| ; scid. flow through nozzle. 
g ere Y is the discharge through a sector subtending one radian at the 
Ss ‘he lial velocity be Ing neglected. the axial velocity uw is uniform 
the « section and is found from Bernoulli's equation, which takes 
29H rgs2 vi. (5.6) 
Surface te hn eiects ove the core are assumed to be imperceptible, sO 
the pi e there is zero. The equation of continuity is 
“a > ~- 
(J (a*—b?). (9.4) 
nsert 1.2), (5.4 | (5.5) into (5.6) and (5.7) and differentiating 
respe the ax ordinate z, we obtain two equations, which 
ite manipuilat n to 
verget May Ski) — #) 
lz R2(2— RY . 
. - (.).3) 
SA(1 R)? 
R3(2— Rs 
da Si 
dz R3(2 Ry 1 
- - (5.9) 
SA (1 R) | 
R3(2 Rys 
Vas s Binnie and Harris (3) that at the throat dw/dz takes the 
mn clven b 
()*(a2— b?) = 
(5.10) 
2h 
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From (5.10) it can be deduced that the denominators of (5.8) and (5,9) a 


required, the values of these quantities can be obtained without difficy! 
by means of the usual method. 


6. Viscous swirling flow through a convergent-divergent nozzle 
It is first necessary to make assumptions concerning the functions ; 


represent wand v. For the former we take 


1 or 
1 a—b (a—b)*} 
where w, is the value of w at the free surface. This expression, whic 
analogous to (3.6), satisfies the conditions 
w O atr a | 


Th w, at? b 
cw 
OQ atr b | 
For v we assume 
= ab? +- (a? +-b?)r—ar??, (5 
b(a—b)? 
which meets the requirements 
v QO at? a \ 


v v, atr b | 
‘ab v 
0 atr b | 


the last condition ensuring zero tangential force at the free surface. F 
bia } the numerical values of (6.1) and (6.3) are displayed in Fig. 
The radial velocity being negligible, the variation of the pressure head 
with radius is shown by Goldstein (5) to be given by the radial equat 


of motion 
cp 


cr r 


with p Oatr b. (6.3 


Writing rs) a (6. 


we then obtain the following results which are wanted later: 
a oC 
% (a) = wy 
Q | rwdr = w, 73 (8a—58); (0 
i = 
a oC 
‘ A «2? : . 
|} ru* dr wy —(16a—118): 
30 


b 


then zero, hence dé/dz and dv,/dz also take the form 0/0 at the throat. | 
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Y ¢ > ‘ . > © »O4 . 
| di UW; (224a3— 483a76 +-368a6?2— 9863): (6.9) 
! 5 1? O0(a—o) 
< , a ' , 
11 2a2(a —5)* log + §(— 12494-42048 
? 6)" | a—o 
NOoZZie 520382 250753 4ad4 69°), (6.10) 
heing the pressure head the wall where r a: 
ye a " o ) 6S 
] . HOdat*(a 5)! log 1 Of 00d -210a°d 
1PO0d4(a—s)* | a—d 
2H0a°d- 12 5a*o” 12a364 42425° 56a0® 158") 4 (6.1 1) 
» 
Wy <3 (6.12) 
oO 
2a—o Tr 
vy (6.13) 
ola 0) 
chee] t was proved that both (6.10) and (6.11) -Oasd>0. 
— T iim ' 
- : 
r ‘ 
| a \ 
5 \ a 
\ 
\ 
\ 
\ 
\ 
o- \ 
5 \ 
L 
\ 
\ 
\ 
\ 
\ 
5 r 0°75 1-0 
a 
Assu distributions in nozzle for b/a }. 
We | » employ bot xial and angular momentum integrals, and 
ixial and angular equations of motion given 
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by Goldstein (5), which for axial symmetry reduce to 


cu cu cp cwu cu cru 
ru PW -— } vi? . dene (6,14 
Ci Cn C v Ci OF A 
; ml 2. : 9 
5 a eo... ct gor cl 9 OU 
and ru ty 1 py = y[7? —-+ vr? —_}. (6.15 
cr Oz or* Cc) oz* 


Here we retain the terms involving the radial velocity uw; but, since th; 
constriction is gradual, we shall follow the approximations of boundary 
layer theory and reject the final ¢? 6z* terms in both expressions. Goldstei 
showed also that the equation of continuity is 


cow 
(ru)+-r 0. (6.16 
cr C2 
Now, from (6.16), 
c c cw cw ow mn 
(raw) w —(ru)+tre ru ru , (6.1; 
a) cr cor Cz Cc} 
cu cw c cu . 
hence ru _7Uur (ruw)-2ru : (6,18 
cr Cz a] Cn 


therefore, if uw, is the value of w at the free surface, 


a a 
Cw ow Pf ” 
| [re Lyrae? dr bu, ws (rae) dr, (6,19 
P or Cz | OA 
b / 
since w = 0 at r =a. Next, we use the theorem that, if f(r,a,b) is 


function of 7, a, and b, 


‘a id di db 
. f(r, a,b) dr : f(r.a,b) dr+f(a,a,b) . f(b, a,b) 7," (6,20 


Cz C2 (le c 


Bo ao ws Y ; » db we 
ru? dr (rw?) dr —bui —, (6.2! 
C2 Cz dz 


again using the condition w= 0 at r =a. From the geometry of tht 


motion it appears that 


(6,22 
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efore with the aid of (6.21) the right-hand side of (6.19) reduces to 
0.14 rw? dr. (6.23) 
‘ g usll 6.20), we have 
a 
( i ; P ( da 
" ] (rp) dr +-ap, ‘ (6.24) 
i G2. J * dz 
Ist e | 0 at 4. Then, on integrating (6.14) and using (6.23) and 


4), we obtain the axial momentum integral 


i eo Ff da Cw _— 
w= di rp dr aps + Ta > (6.25) 
( C2 az cl 
: r=a 
being zero at ? b. If we consider two cross-sections of the nozzle 


small distance dz apart, we find from first principles that the physical 
rpretation of (6.25) is as follows: the increase in the flow of axial momen- 


jual to minus the increase of axial pressure over the annular section 


~ { 
s the axial component of the wall pressure minus the viscous force due 
The angular momentum integral will now be derived by a similar process. 
t a) 

( ¢ eo. Ga o.. 08 . 

? ) u (rv) —7°V | 7-4 + ruv, (6.26) 
( cor C2 cr 

. Cov , ? c 9 Cc 9 > = 

nee “1 { L yuy (r-uv) + (r-vw), (6.27) 
( cr C2 
elore 
a 

‘ € ) uC 9 ; c » » , 

{ rw - uv) di —b*u, vy (r-vw) dr, (6.28) 
5 : 
is zero at the wall. From (6.20) and (6.22) we have 
¢ » c 9 9 > ) 
r-vw d (r-vw) dr bu, V1; (6.29) 
C2 
gain using the condition that v 0 at r =a. Thus the right-hand side 
6 IR) 3 oe, 
- IS SIMpl\ 
a 


(6.30) 
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Now 
ae PL i cv > ev i ia ,| ov 
| gi tr - v\ dr f* — —17v a ‘ (6.31 
5 iia cr or or 
. b r=a 
b : 
because at r= a, v 0; and at r=), dév/ér v/r in conformity with 


(6.4). Then, on integrating (6.15), we arrive with the help of (6.30) and 
(6.31) at the angular momentum integral 


a 
> 


row dr va 





ov 
. : (6,3? 
cr 

r=-a 


Again, this may be expressed in words as: the increase in the flow of angular 


b 


momentum is equal to minus the viscous couple exerted by the wall, 
The term omitted in (6.15) evidently represents the difference between the 
viscous couples exerted on the annular ends of the element by the neigh- 
bouring layers. 

The formulae (6.7)—(6.13) are now inserted in the momentum integrals 
First, w, is eliminated from (6.8), (6.9), and (6.12) by means of (6.7), so 
that the variables, with which we shall be concerned, are 6, a, and ¢,. 
To obtain results in non-dimensional form, we again use (5.4) and (5.4), 
and add 

i= =. (6.33 
Q 
For the angular momentum integral, (6.9) and (6.13) are required, fro 
which, after performing the necessary differentiation, we obtain without 
further approximation 
dd a dv 


B ' A da 


_CL, 6.3 
dz v, dz dz | 


woe 


where A, B, and C, which are functions of R only, are given by 


1792— 5824 R+ 6695 R2— 3262 R?-+ 566 R! 


(S—5.R)(1— R)(224—483 R-+ 368 R2— 98 R3)’ 


B 952— 3648 R-+-4721 R?— 2548 R3+-490R! (6.3 
(S—5R)(1— R)(224—483 R-+ 368 R2— 98 R3)’ 


C 36(8—S5HRU2— RR) 6.37 
- 5 (0.01 
R(224— 483 R+ 368 R?— 98 R3) 

The axial momentum integral involves (6.8) and (6.10)—(6.12), and fron 

them after very considerable labour we find 
Pee dd a dv, , ., da j 9 
(EK —H)-" —F 7 = (G—DK) + KL, (6.38 

: r 


Ge VU, de 








in wl 


Th 
dv 





Thus we have two simultaneous equations (6.34) and (6.35) for dé/dz and 
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n which 


D 32R(4 a) (6.39) 
d(8 5k)? 
. 2(64—120R-+-55 he) (6.40) 
5(8—5 Rk)? 
8—5R 
F ‘__ 60(1— R)*log(1— R)+ R(— 60+ 210R— 
1440 R2(1— R)2 


260 R2+- 125. R3—12R*4 42 R®5—56R8+15R’)}, (6.41) 
G 85K ¢__ 601 R)(2—R)log(1— R)-+ R(—120+480R— 
1440 R2(1— R)?* 
730 R?2+-510R3— 181 R4+- 68 R°—44R6+13R")}, (6.42) 
H ; 1s 60(1— R)4(2— R)log(1— R)+ R(120—480R- 
1440 R3(1— R) 
730. R2— 510 R3+ 149 R?2+ 36 R5— 84 R64 58.R7—15R8)}. (6.43) 


dz. The results of solving them are 


] 


dd ‘DA -+AF—G}(da/dz)+{CF—K}L , 
— —— : (6.44) 
adv, ((AE+ BD)kK—(AH+ BG)}(da/dz)—{(B—CE)K-+CH}L. 
v, dz EK —(BF--H) 
(6.45) 


hey should be compared with (5.8) and (5.9), which are the corresponding 
results for inviscid flow. Again, dé/dz must take the form 0/0 at the 


effective throat, the position of which is therefore given by the two 


conditions EK (BF | H) 0. (6.46) 
_ _, da _ : a 
{‘DA+AF—G! | {CF— KIL 0. (6.47) 

dz 
On inserting (6.46) into the numerators of (6.44) and (6.45), we find the 
latter to be equal to the former multiplied by B, without reference to 


6.35)-(6.37) and (6.39)-(6.43). Hence dv,/dz also takes the form 0/0 at 
the effective throat. The result (6.44) is so complicated that it is im- 
practicable to derive an expression for the throat value of d3/dz by the 
usual method of differentiation. 

It will now be proved that a long wave of small amplitude can maintain 
itself stationary on the streaming and revolving surface of the liquid at 
the effective throat. As shown in Fig. 5, we suppose a stationary wave to 


exist having an elevation « at the section XY so that the surface there is 
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at radius b—e, and we denote by w’ and v’ the small additional axial and 
tangential velocities at Y due to the wave. In conformity with (6.1) and 
(6.3), the distributions of axial and tangential velocities across XY are 
taken to be 


nf 2(a—r) (a—r)? |) : 
i] (u w _ — — . 6.48) 
d 1 Na (b—e) fa (b—e)}?| (6.48 
v i | —a(b—e)?+ {a?+-(b—e)*}r ar}. (6.49) 


(b €)\a (b—e)}? 


W444 Vt 


i phn oll 
rs - - 2 
€ b-e * ¥ b 


se y t= 


Fic. 5. Wave on free surface in nozzle. 











A number of integrals across XY are required, and they can be obtained 
from (6.7), (6.8), (6.9), and (6.11) on replacing b by b—e, 5 by 6-+-e, w, by 
w,+w’,and v, by v,+v’. Thus with the substitutions (6.48) and (6.49) we 
obtain 


a 


‘ ’ d+e ~— — 
Q |} rwdr (w,+w ( 13 84-36 }-¢)t, (6.50) 
b—e = 
a 24, 16a—11(8+< - 
ru dr — ¢)* = ~ ) 73? (6.51) 
cs 5 (5+ e)(8a—56-+e)} 
- ‘ ] @ | 12240" 483a?(d-++-e)- 368a6(6-+- 2e) —988?(6 + 3e)! 
r-vw di —(U v - eT : oe = 
ae = {sa—5(d-+«€)Ha—(d+e)} 
(6.52) 
a (vy, a’ )2 : . a 
| rp dr - u ———. | 60a*{a — (6-4 <)#log! ea 
r 1206°(6-+-4e)fa—(d+.e)!" la—(d+e)} 
60a7(5-+-€) + 210a®5(5 + 2e) —260a°d?(5+- 3e)- 
125a483(6 +-4e) — 12a364(8 + 5e) + 42a75°(5+ 6e) 
56a6°(6-+-7e)+1567(6+-8e)]. (6.53) 


In (6.51) and (6.52) w,-+-w’ has been eliminated by means of (6.50). The 
axial momentum equation (6.25) is now 


a a 
C F ? C ‘ rr 
rw? dr rp dr. (6.54) 
oz az . 
b—e b-e 


We therefore differentiate (6.51) and (6.53), remembering that v’ and « 
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and 6 are constants. After protracted calcula- 


nd products of small quantities are again ignored, 


Thies 
24Q0°E Le 24av7H de 24a*v,F dv’ (6.55) 
‘ 37) 
}R2(8—S5R) dz R2(8—5R)dz R(8—5R) dz 
‘ de a dy’ ° 
luces to EA —H) ? (6.56) 
dz v, dz 














1 
1:0 
‘tive throat of nozzle: viscous, - - - - inviscid. 
ir momentum equation (6.32) in this instance is 
a 
r2vu dr 0. (6.57) 
6.52) has been differentiated, yields 
dv v de 
B-! (0). (6.58) 
dz a dz 
6.56) we arrive at the condition 
EA —(BF-++-H) = 0, (6.59) 
ent with (6.46). There is, of course, no result corresponding to 
ch involves 1 | da dz. 
tions that exist at the effective throat may be stated in another 
Equations (6.46) and (6.7) are expressions for @Q?/(vja*) and 
hence by division the interesting ratio v,/w, is given by 
R -~D E \3 - 
(3d ott) ~ . (6.00) 
Te 12 |\BF--H/ 


has been plotted ona base of R in Fig. 6, where it is shown by 
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the full line. The combination of (5.1) and (5.10) yields the corresponding 
result for frictionless flow which is 


(6.61) 

















0 2:5 5 75 P 10 


Fic. 7. Profiles of nozzle and stream: viscous, - - - - inviscid. 


7. Numerical example of swirling flow through a convergent- 
divergent nozzle 
The nozzle is taken to be of the symmetrical parabolic shape defined by 


a 0-4 0-082(I : ). (7.1) 
10 


It extends from z = 0 to z = 10; and at the throat a = 0-2. This shape 
accords with the requirement that the constriction should be gradual: 
it is shown to very unequal scales in Fig. 7. The conditions at the effectiv 
throat are assumed to be z 5:3, R = 0-5, and these lead to a negative 
value of d3/dz which is wanted in order that comparison may be made 
with the free flow of inviscid liquid. These conditions were not chosen 
without some preliminary trials. The value of dd/dz is obtained from a 
quadratic equation analogous to (3.15), both roots of which may be positive 
if the viscous forces are made excessive by the selection of too small a value 
of R. For the above values of z and R it follows from (6.46) and (6.47) that 
K 0:255327, L 0-00821695. If we suppose (7.1) to be expressed in 
inches and take the liquid to be an oil for which v = 57-6 10-3 in.?/sec., 
we find from (6.33) that Q 1-407026 in.?/see. 

The solution was very kindly obtained by Dr. J. C. P. Miller and 
Dr. D. J. Wheeler. The two main difficulties were due to the 0/0 form 
of the equation for dd/dz at the effective throat. The first difficulty arose 
in evaluating the numerical values of the coefficients of the quadratic 
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equation mentioned above. These depended on derivatives of very compli- 
cated expressions which it was not feasible to differentiate algebraically in 
full. Instead, numerical derivatives were computed for each intermediate 


function used in the evaluation of the expression; that is, algebraic substi 


tution was not carried out completely before numerical evaluation was 
started. As a simple example, if the derivative of 
Yy (23 +-32x)log x 


was needed, the expression 


LY 


(3a 3)log x (a7?--3) 


would never have been obtained: instead, numerical values would have 


been found fol 


, ets 4 , I 
r? + 32 { oxre°+ 3, v log x, v : 
x 
ind then y uv’+-u'v would have been evaluated. Though the variation 


in method is trivial in this example, the alternation of computational and 
ugebraic processes in this way made a considerable difference to the ease 


with which computations were carried out. 


The second difficulty arose in the evaluation of the solution to 


d8 /dz ND, say, 


in the neighbourhood of the effective throat, where NV D = 0, so that 
t was very difficult to obtain an accurate value for the quotient. This was 
overcome by solving two extra differential equations in the neighbourhood 
of the throat, that is dN dz N’, dD/dz D’, where N’ and D’ were 
computed as outlined above, and by using these equations to give N and 
D to the desired accuracy. 

The calculations were performed by EDSAC. The machine provided 
values of R and K at small intervals of z in both directions along the nozzle, 
and by means of (5.4), (5.5), and (6.7) the values of 5, v,, and w, were 


obtained that are shown in Figs. 7 and 8 by full lines. For comparison, the 
broken lines in these figures give the corresponding quantities for inviscid 
flow on the assumption that the discharge and the core radius at the centre 


of the hozzZile are the same as before. 
In this example the effects of viscosity are large. With inviscid flow the 
thickness 8 so rapidly diminishes from inlet to outlet that the core radius b 


increases throughout the nozzle. In viscous flow the reduction in 6 is much 
less marked, and near the outlet, over the range 8-8 < z < 10, the opposing 
lorces are sufficiently great to cause a very slight increase in 6, which is 


hardly noticeable in the diagram. The changes in thickness are so small 
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that in the convergent part the core radius decreases to a minimum neat 
4:4. Fig. 8 shows how the inviscid values of the tangential surface 
velocity v, vary through the nozzle inversely with the core radius. Tp 


viscous flow a larger value of 7, must be imposed at inlet; the reduction js 











Kr <a ompone nt velocities at fre« surface mn nozzle: VISCOUS, - - INVIS 


now sharper, and this component of velocity is almost extinguished when 
the outlet is reached. The streaming velocity w, in inviscid flow, which is 
uniform over the cross-section, increases throughout, but in viscous flow it 


reaches a maximum near z 5-2 and thereafter diminishes. 
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{ NOTE ON THE OSEEN APPROXIMATION FOR A 
PARABOLOID IN A UNIFORM STREAM PARALLEL 
TO ITS AXIS 


B j | \\ | LK | NN yN (Unive rsity Colle qe o} Swansea) 


Received 24 March 1955 


SUMMARY 
the Oseen equations have a simple solution when 
loid of elliptic section in a stream of viscous, 


tv has uniform velocity in the direction of the axis. 


yn is ! the case in which the paraboloid degenerates into a 
rab vding edge, and it is suggested that the solution of the 

seer tl : vy throw some light on the corresponding three- 
lar 


1. Introduction 


ARRIT | Lewis (1) first gave the solution of the Oseen equations for 
half-plane QO, 2 Qin a uniform stream. At about the same time 
Stewarts working independently, produced the same solution by a 


imple method which also gave the solution in the case of a 

paraboli inder. Stewartson did not publish this result. A recent paper 
Kaplun (2) contains a solution similar to Stewartson’s. 

For the parabolic cylinder the solution of the Oseen equations can be 


In the form 


/ 0) u cp eX (1) 
’ 2 ( oh CZ 
where mee ei OX, EX _ op OX k = U/2 
i! “ya ad Ca 

nd wu 0 at infinity and wa Uw 0 on the cylinder. 

Stewartson expresses these equations in terms of parabolic cylinder 
oordinates €, » defined by x = €*—y?, z = 2&y in which the cylinder is 
the coordinate line onstant. He then seeks, and finds, functions 4, x 
which are functions of nly 

The question natw lly arises: Can this method be extended from the 
parabolic cylinder in two dimensions to the paraboloid in three dimensions 


by the use of paraboloidal coordinates ? It isshown below that this extension 


Quart. Journ. Mech. and Applied Math., Vol. VIII, Pt. 4 (1955)] 
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is possible for paraboloids of elliptic section with vertices facing into the 


uniform stream. 


2. The equations 


Paraboloidal coordinates (€,,¢) can be defined as the roots A of the 


equation 


2) > 
y- no 
* + A(x -r) (? 
a—A ° b—A 
when x, y, 2 are supposed given, a > b > Oandé >a,a >yn>b,b>l 
Then the surface € = ¢, is the paraboloid of elliptic section 
y" ~2 
4(x— Cy) 


and by asuitable choice of coordinates in the first place any such paraboloid, 
except that with a circular section, can be expressed as € = C, if a, b and é 
are correctly chosen. 


The coordinates x, y, 2 are given in terms of €, », ¢ by 


24 (€ b)(y b)(b—Z). 


[f the flow has velocity U in the direction of the a-axis of infinity, the Oseen 
equations are 


eu | Op Cu cu atu 
| Fae Te. 
COX pcx ox” cy~ CZ° 
, OV le¢ p Cn en en » 
l v( ) (5 
Ox p cy on oy* 6s" 
ew l ep ctw ctw *w' 
Ox p Cz cy= ay" 
ou al) ou 
0 +. - a 
OX cy Cz 
with u —U,v=w=0o0n ¢, and uw, v, w > 0 at infinity. 


In the general three-dimensional problem it is not possible to write the 
solution of the equations (5) in a form corresponding to equations (1), but 
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n this case a solution is to be sought of the form 
Cd ] ox 
C2 2k O02 x 
Ch l CX 
oy 2k ey : 
; J (6) 
Ch | ox 
C2 2k ez’ 
od 
pl ; 
Cx 
ond o*d o*h P= 
en \ The Fin fe (7) 
Cx” oy” Oy Aa 
C= 24, Czy A Cc , 
v2) S| — 24. om SEE, (8) 
Oox*  oy* d2* Ox 
vhere 2v and in paraboloidal coordinates 
—— 
bs &€—a)(E—bD c f , 7. 
\ ¢ a)(§—b)}* — 
(é Che ? a o€ 
(a (7 byt? ¢ C 
{(a—n)(n—b)}* —| + 
(7 Sis—3) ont on 
(a—C)(b—Z)¥ ¢ 7 soa 
oa 1 {(a—f(b—f)}+# — 
Y—Ghe—6) OG as 
ay(E—b) « (a—n)(y b) « (a—l)(b—L) ¢ 
(E—n) CE (E—H)(H—S) On  (y—LYE—f) oC 


The existence of such a set of solutions implies that the vorticity always 

es in a plane perpendicular to the axis of the paraboloid. Moreover, the 
rm of equations (7) and (8) in paraboloidal coordinates makes it possible 
find solutions of these equations which are functions of ¢ only. 


In this case (7) and (8) become 








d | 
r—Ll)(b—C Ad = 0, (9) 
d¢| ; : dl 
d I dx O16 y a 1, 4x 
I S(a—C)\(b—2 2ks(a—Cy(t pA (10 
dt | dt 1 ( S)\Y G)5 dt ) 
Solutions of these equations may be written as 
¢ = Alogia/2+b/2—¢+{(a—L)(b—Z)}], (11) 
x e2ks 
Y B i ds, (12) 


‘(a—s)(b—s)\ 


where 





> constants. 
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Equations (6) can now be written as 


(a 6)(6 6) “(6+ | ) . 


(E—CL)(n—C) dC 2h* 
((a—L)(b—L)}¥ ( , , Bers) 
—— {Ao x: 
(E—C)(y—S) | 2k | 
‘ (b- C)i(gé—a)(a n)(a -C)rt “(¢ | 
(a—b)'(E—CL)(n—Z) dt \" 2p) 
(13 
{(é—a)(a n)(b C)M { r Be 2ks) 
(a—b)'(E—L)(n—L) 2k f’ 
(a—l){(E—b)(n—b)(b—L)}' d 
w 7 a P+ X 
(a—b)'(E—L)(n—f) dc 2k” 
(é b)(y b)\(a—L)4 { r Be au 
(g—b)"E—L)(m—) ( ° 2k J 
22 
The conditions » = w = Ofor ( = C,are satisfied by taking A = 
2h 
and the condition w U is satisfied by taking 
I fe a 
(Co) > Ss fa 14) 
ai ((a—s)(b—s)}* | | 
It is clear that wu, v, w > 0 as €-+ —o which corresponds to moving off 


to infinite distance from the solid, and it is not difficult to show that th 
vortex lines are the plane ellipses 


C constant, x = constant. 


Equation (3) cannot represent a paraboloid of circular section unless 


a — 6b, when the coordinates become degenerate and the method fails 
However, a result due to Howarth may be employed in this case. The 
coordinates €, 7 are introduced, where 


: i 9 e ) i Ly 9 
a "Seca | y*-+-2* 4&“n°. 


A solution of the type given in equations (6) can then be found with ¢ 


and y functions of 7 only. The equations corresponding to (9) and (10) ar 


the equations d { dd 
(» 0 (19 
dy dy) 
] ly ly 
i (aa) = ae 7 
dy \ dy ay 


. , ¢ ‘ia 
with solutions d = A log », eae 5 | ds. (li 
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The boundary conditions can then be satisfied by taking 2/A Be-?k 


here the paraboloid of revolution is the surface 7 No- 


The solution of the Oseen equations specified in equations (11), (12), (13) 


nd (14), will in particular apply to the degenerate paraboloid € = (, = 6, 
which is a flat plate with a parabolic leading edge. In dealing with this case 
t is convenient to take 6 = ¢, = 0, when the parabolic plate is z = 0, 
/ 
B—U| | — ds Ueka Vek 
I > = - 9 
| a—s); K, (ka) k,( R/4) 
here R is a Reynolds number defined as 
2aU 
R tak, 
* 
being the semi-latus rectum of the parabola. 
lhe fluid velocity u given by the approximation may then be written 
u = U+u,-+u4,, 
here U,0,0), u, y,0,0), and uy, is a vector of magnitude 
2aUe-*A l—e-? ae ; 
directed in the ¢-direction—.e. along the tangent 
RK,(.R/4 (7 
) the curve ; constant constant in the direction of ¢ increasing. 


e vector u, influences only the 2-component of the flow and is a pure 
oundary layer term. As # increases it is significant in a smaller and 


! 
y 


smaller region round the plate. The secondary flow—the three-dimensional 


fect induced by the presence of the plate—is given by the term u,. The 
gnitude of this term is zero on the plate and at infinity and decreases 
el Ri as R-> x. The direction of u, is always away from the plate and 
the limit, as the plate is approached, is normal to the plate. 


This secondary velocity u, may be considered as due to a source distri- 


Ue-F4, C(a—Z): , 
ion of density = per unit volume; and this density 
K,(R/4 C)\(n—C)} 
ilwa ILIVEé 
\ «> 
Again this case it is a simple matter to show that on the plate 
Cu B ov 
| : | 0 
Cz 0 C7) C2] ~-9 
t igh the flow is strictly three-dimensional, the skin friction is 


the direction of the flow at infinity. Moreover, when ¢ 0 equations 
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(4) can be used to show that &y \(4ax—y?) so that 

ou 2B Ueki 

(: Z 0 (tax y*)? ae Kol R 4)’ (18 
where 7 (a—y?/4a) distance from edge of plate measured parallel t; 


the x-axis. 
It is interesting to compare this result with the values of (éu éz 
obtained from other methods of approximation. If the boundary-laye 


equations are taken in the form 


cou ou ou o-u 
u i uw V 
> 
Cx oy C2 C2* 
‘ A 
Cv ol ol ov 
u a w 1 = (19 
On cy CZ C2" 
cu. cv, ow 
U t | 0, 
Ox oy Cz 


Moore (3) has shown that there will be a solution with 7 0 and on the 


s 0-332U R! 
3] 0 


9) sh 
The ratio of this value to the value given by equation (18) is independent 


plate 


aa 


of position on the plate. 
However, the simplest approximation to the skin friction on the plate is 
obtained by taking Rayleigh’s value for the skin friction (4), starting the 


(™ UR} - 
‘ (21 

a 34444 -04 

z),-9  2atatxt 


\ z 


solution from the leading edge of the plate. The value of (éu/éz)._, obtained 


1S 


The ratio of this value to the value obtained from the Oseen equations 
in equation (18) is independent of position on the plate and is a functio 
only of the Reynolds number R. The way the ratio varies with F is show! 
in Fig. 1. Beginning at zero for R = 0 it is greater than 0-9 when R = 4 
and rapidly approaches its limiting value of unity. Thus, except for ver} 
small values of the Reynolds number, the Oseen equations provide almost 
no information about the skin friction that cannot be obtained far more 
easily from this simple Rayleigh approximation. Nor is this result altogether 
surprising, for the equations which determine the Rayleigh approximation 
are obtained by forming boundary-layer equations like equations (19) from 
the Oseen equations instead of the full equations of motion. What is 0! 
interest is the fact that the omission of the é2/éy? terms from these equations 
should have so little effect at large distances from the axis of the parabola 
for it would be expected that far from the axis, where the direction of the 
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‘tangent to the leading edge approaches the direction of the axis, ‘edge 
fects’ would become important and the terms in @?/éy? would not be 


egligible in comparison with the ¢?/¢éz? terms. 








_ 
i, 
a4, 
lt 
4 
+ — 4 4 4 4 1 -® i 
4 12 16 20 24 28 32 
Reynolds number R&R ——> 
Ky l. Ratio of skin iction from Rayleigh approximation to skin 


from Oseen equations. 


lhe fact that edge effects are of so little significance in this approximation 
ggests that equations (19) may give a valid approximation to the skin 
tion over the whole plate. Equations (19) ignore edge effects, but as 
ey bear the same relationship to the full equations of motion as the 
lations determining the Rayleigh approximation do to the Oseen equa- 
ns, it seems reasonable to suppose that edge effects will have very little 
fect on skin friction for high Reynolds numbers. If this is the case, 


») 


juation (20) should give a much better approximation to the skin friction 


in equation (18) for high Reynolas numbers, for the limitations of the 


seen equ itions are well known. 
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THE INITIAL VALLE PROBLEM FOR THE WAVE 
EQUATION IN THE DISTRIBUTIONS OF 
SCIIW ARTZ 


By J. OKRKEFFE (University College, Galway) 
{Received 17 February 1955 


SUMMARY 
The classic formulae for the solution of initial value problems for the w 


equation are found in the ease of distributions by a Fourier transform method. 


1. Summary of classic cases and results 
THE following are the classic formulae for the solution of the homogeneous 


wave equation 


1 ctu 
Liu) Au = — 0, A 
C2 
with initial conditions 
c 
u(r, O) U(X), u(x, O) u,(v), 
ct 
k ‘ 


. ° oe Cg . . 
where A is the Laplacian operator S —, in k-dimensional space and 2 is 
m— OF 


i] 
the point (2,2 ,...,2,): 
(i) the case / | (equation of vibrating strings), d’ Alembert’s solution 


a et 


> 


l ee 
u(a. t) YW ag (a ct) u(x ct)|4 : | wi (&) dé: 
aC 


J ~ ct 
(ii) the case 4 3 (equation of spherical waves), Poisson's solution, 
d ‘ 
u(a, ¥, 2, t) ig Lt Malo) | -tM(u,), 2 
( 


where J/,(w) denotes the average value of wu over the surface of the sphere 
of radius r with centre at (a, y, 2): 
(iii) deduced from this for the case k 2 (equation of cylindrical waves 


the Poisson—Parseval solution, 





1 [d 9 
u(x, y,t) — qy tect (to) + phg( Uy) | > d 
rf ul(€, ; 
where [,.(w) | | — ae er dédy. 
E244 n)t<r? Vi . C é W 3 


(Quart. Journ. Mech. and Applied Math., Vol. VIII, Pt. 4 (1955)] 
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e original derivation of Poisson’s solution consisted in verifying directly 
it the formula (2) satisfies every required condition; a theorem of 
Imgren then proves that the solution is unique (see 1, 47). 
iv) the case « l,/ 2m-+-3, Tedone (2) found a formula equivalent 
that found later by Hadamard (1, 244) as a particular case of his ‘finite 
rt’ solution of the general non-parabolic partial differential equation of 
second orde 
k—3 
; vt) “i th 3M)(U,) | d | 3M(u,) |, 
dt\2td 2tdt 
ere s the area of the unit hypersphere in k-dimensional space and 


uw) denotes the average value of u over the hypersphere of radius 7 with 


ntre at Hadamard’s method also furnishes a direct proof of the 


<istence of Poisson’s solution. The same solution was found in the 


ternative torm 


| (e p2) hh My) dr (4) 


(a1) 


0 
Courant (3, 387) by a series of formal operations on multiple Fourier 
ransforms, the resulting formula being verified a posteriori as satisfying all 


required co..ditions. From this solution he then deduced by Hada- 
rd’s method of descent 
» 


the same formula (4) in the case k 2m 


Courant’s method was subsequently made rigorous by Cooper (4) by 


he use of spherical Abelian means. Later, Bureau (5) found formulae 


juivalent to (4) for the odd and even cases separately by generalizing the 


lefinition of the Fourier transforms occurring in Courant’s and Cooper’s 


ethods to include the finite and logarithmic parts of divergent integrals. 
The classic formulae for the solution of the non-homogeneous wave- 
equation 
Liu) q(x, t) (B) 
th 


h homogeneous initial conditions are 


uiorft be ( [ glé, r) dédr, (5) 


where the integral is evaluated over the triangle bounded by the straight 


ies 


c(t T) (a &) 0, T Q: 
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(vii) in the case k = 2, 


whet, t) =e © | | | . g(€, 9.7) - dédndr, (6) 


27 c?(t—7)?— (a—€&)?—(y—7)?! 
where the integral is evaluated over the region bounded by the cone 
c?(t—r)?—(x—£)?—(y—n)*? = 0 


and the plane 7 = 0; 


(viii) in the case / = 3, the retarded potential, 
1 (ff g(é,n,6,t—r/e 
u(x, y,2z,t) = | | | HE 76 9) dédndt, (7) 
4a J J. r . . 
r=ct 
where r? = (~a—&)?+(y—n)?4+(z—L)? 
(ix) in the cases k > 3, taking c = 1, 
wine @) ] (2) 2 (t2 — 2) MM g(r,t)] dr (8) 
(k—2)!\et ‘ - ; 
0 
ke 
where ry? = } (2;--&;)? 
i=1 


(3, 402).- 


2. The Fourier transform method 

Courant’s method is essentially a method for calculating the elementary 
solution of the partial differential equation (A) and then using what 
Schwartz (6 (ii), 67) has since shown to be a general property of convolu- 
tions of distributions. This property ensures that any partial differential 
equation can be solved once the elementary solution is known. Gates (7 
has applied it to the case of the initial value problem for the equation of 
heat conduction in one dimension. 

The elementary solution has already been given by Garnir (9 and 8, 278 
for the cases k = 1, 2, 3, and has been given by Schwartz (cf. 6 (ii), 146) for 
all k in the form 

i =! ie 

——, F% Pf, oe = of—s%, 

4a“ o* 

where ¥ denotes the inverse Fourier transform in the sense of the theory 
of distributions with respect to both space and time coordinates, and Pf . 
denotes the distribution corresponding to the operation of taking Hada- 
mard’s ‘finite part’ of the integrals of products of o~? and test-funct ons. 
Garnir finds his elementary solutions by a Laplace transform method some- 
what similar to the Fourier transform method used in this paper. Ga:es 
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merely verifies that a certain distribution is an elementary solution. In 


Bureau's method for the classical case of functions the elementary solution 


s the same as that of Hadamard, but this is not necessarily the same as 


the elementary solution in distributions (cf. 6 (i), 133). 


Following Courant’s method, in this paper the elementary solutions are 
ilculated directly by considering Fourier transforms with respect to space- 
ordinates only and in all cases the formulae 
sin 2zcpt 
) ? P 


ke 
> 
«77 p fan ‘ 


1 


found, where .Y denotes the Fourier transform with respect to space- 


ordinates only. This expression is easier to evaluate than Schwartz’s. 

The case of distributions subsumes all the classic cases so far studied, and 
he advantage of dealing with them instead of with the classical functions 
s that all difficulties with regard to convergence and differentiability of 
ntegrals are eliminated. Thus the formal operations of Courant are imme- 
liately justified and the cases of odd and even numbers of variables become 
similar. As compared with the method of Bureau the method of distribu- 
ons has the advantage that it leads directly to the final expressions in terms 
f convergent integrals without having to consider finite and logarithmic 


parts. 


3. The theory of distributions 


The theory of distributions (6) is essentially a theory of weak convergence 


nd for our purpose may be summarized as follows: 

\ distribution is a linear functional on the space (7) of infinitely differen- 
tiable functions of compact support. Test-functions (‘fonctions indéfini- 
ment dérivables a décroissance rapide’ (6 (ii), 89) are defined by the 
onditions that 


i) every test-function ¢(x) is infinitely differentiable at all points x in 
dimensional space 

OX+B+... +4 
ll) lin 14 X5 ae - 5 P(X) 0 


lor all positive integers /, m,..., DP, O, By.oy pe 
A parameter function of slow growth with respect to « (‘fonction tem- 
perée’ (6 (ii), 93)), i.e. any function of the form 
x |*)**f (a, ax)}, f(x,a)| < K, 


( x4 ¢ NG 00 CQ] 


where all the indices are positive integers or zero and k > 0, is said to 
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converge to the limit /(x) in the space of distributions of slow growth as y 
tends to X45 if 


lim ( | f(a, ) —Ua) |\b(ar) dx 0 


in the Lebesgue sense for all test-functions 4(2), the integral being taken 
over the support of (x). In particular, the Fourier transform of a bounded 
function f(a) 
. k 
F f(x) = | e®*"*$f(x) dz, a. => 2, €,, 
. 1=1 


converges in the sense of the theory of distributions to the limit F(é) if 
| [F f(x)—F(€)|d(é) dé = 0 


for all test-functions. 

The principal properties of distributions used in this paper are the 
following: 

(i) A distribution of slow growth has an infinite number of partial 
derivatives, where the operation of differentiation is defined for distributions 
as in (6 (i), Ch. 2) and, in the cases where the distribution is a function 
differentiable in the ordinary sense, the derivative as defined for distri- 
butions is identical with the ordinary derivative. In particular, an integral 
containing a parameter may be differentiated with respect to that parameter 
if it can be shown to converge in the sense of the theory of distributions. 
This property will ensure the validity of all the operations performed o1 
Fourier transforms in Courant’s method. 

(ii) The Dirae delta function 


{ 0, le i 
O(2x) 


| d(x) da :; 
x. £ 0 4 
and its derivatives are ‘measures’ and may be integrated and differentiated 
without restriction when associated with any continuous function of con- 
pact support (6 (i), 11). 

(iii) If the Fourier transform of a function converges in the sense of the 
theory of distributions, the original function is the unique inverse Fourier 
transform of that transform (6 (ii), 107). 

(iv) If Lis the distribution corresponding to a given differential operator 
and the distribution F satisfies the differential equation LE = 4, where 
5 is the Dirac measure in the space of all the variables concerned, the 
distribution EF * B, defined in the case where E and B are integrable 
functions as 


Ex B= | dr | dé, [ dé... | dé, Ele—é, t—7)BlE,7) 
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nd known as the convolution of E and B, satisfies the differential equation 
LT B, where B is any given distribution (6 (ii), 6 and 67). When B 
snot an integrable function the definition of FE * B can be suitably general- 
ized but this will not be needed here. 

v) Garnir (10, 90) proves that where the solution of a partial differential 
equation in the sense of the theory of distributions satisfies the requirements 
for a solution in the ordinary sense it is identical with the solution in the 


rdinar\ ense. 


4, The elementary solutions 


\n elementary solution of either of the differential equations (A) or (B) 
1 the theory of distributions is any distribution EF which satisfies the 


relation ’ 
LE Oo, 
where the Dirac measure for the space of all the variables occurring 
the differential equation. 
Assume (following the method of Schwartz (6 (ii), 142) but taking 
Fourier transforms with respect to space coordinates only) that 
E(x, t) FEE, t), 


where * denotes the inverse Fourier transform with respect to space 


ordinates. By differentiation under the integral sign we have 
| l &é - - 
LE j br 9%¢ 5 sment 5(x)d(t), p* == > £7, 

C* a iz 


] 


ind, taking the Fourier transform of this equation, we have 


> L 


$7770? 076 c76(t). 


This differential equation can be solved by the method of variation of 


parameters, giving the elementary solution 


, : sin 2zcpt 
E(x,t) = F&(E,t),  Elé,t) ie, 


7p 


where 
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When k = 1, 
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7 _ 2rcét i ¢ SI 27rcé 
B(x,t) = —cH(t)F¥ 2 — — chit) | e-2nige meet ge 
ome . 2ré . 
= 
5 H(ct vr) 


by contour integration. 
When k = 2, 


. 


H(t)F E(é,t) = —cH (t) | pdp o-27irp cos 951 27Cpt 16 





i 2arp 
0 0 
—c?H(t) | sin 2zcpt J,(27rp) dp 
0 
c H(ct—r) 


2a ./(c%t2?—r?)’ 


by Weber’s formula (11, 405). 
When k = 3, 


I. 7 





z ' , , sin 27¢ 
H(t)F &(é,t) 2ncH (t) | p2 dp | e~2rirp cos 90 a do 
a7Tp 
: ; : 
H(t) (. 
Cc (2) sin 27cpt sin 27rp dp 
2rp . 
0 
H(t) ¢ 
mil (t) | | cos 27(ct —r)p—cos 27(ct+-r)p]| dp. 


27p . 
0 


The integral is divergent in the classical sense but may be calculated in 
distributions as follows. Consider the integral 


R 
. = a 

lsin wp| sinaR 
cos ap dp = p (a. == 0), 
“ 0 x 

0 


and also the integral 


x 


I(R) = | $(a)sinaR de, 


x 
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where (a) 1s any test-function. Now 


4 Lr /CUOS R : | , , 
1(R) : d'(a)cos aR da 
R R 


. 


nd therefore lim /(R 0, since 


vsinaR da ( d'(a)\da < kK, 


L 


so that. in the sense of the theory of distributions, 


( 0, ifet—r ~9, 


Hit | |}cos 2Prict } cos 2a(ct rj)p dp aia 
pl ds \—oo, ifct—r 0. 


ck 


H(t) F &é(é,t) d(ct—r) 
Dany 
= 
ind K IS @ASLI\ found to pe 3 by considering the inverse relation 
oi sin 27c pt 
“ (cl—r) = C€ ; 0. 
oT 27rp 
\\ ) : heretore 
' c o(ct—r) 
k j 
47 ? 
san eleme nt ry solution 
When f 4 
2a)?" * sin 2rcpt 
t ( bk re 
, Dany ik —2 | = pS iy »)(27rp) dp 
- . =7 p 
0 


»)(27rp) dp. 


sin 2rc pt pik Wri 


0 


Again the integral is divergent in the ordinary sense, but may be shown 
to have a meaning in the theory of distributions as follows: we have 


1 


D(- k-2 * 


és >) o P 
J (2 an, awe NIG (1—u?)i*—eos 2arpu du 
irs). 

0 
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by Lommel’s formula (11, 48). Now, before proceeding to the limit in the 
formal expression 


R 1 
k—2) . . 
i). | p* 2sin 27c pt dp (1—wu?)4—Deos 2arr pu du 


. 


0 


—lim 


lim naa 





consider instead of the double integral the equivalent double integral 
1 R 
I(R) | (1—w?)i*-3) du | p* 2 3| sin 27(ct—ru)p+sin 2z(ct4 ru)p| dp. 
0 0 
The inner integral is of the form 
. R 


COS xp R Pp 


| p” sin xp dp | — p? p”-l cos ap dp 


e | xX lo x 
0 0 


vate R 4 I 
p COS xp 4 P p , un xp 


p 


JP p” 2 sin xp dp, 


and so on. In the theory of distributions 
lim p’sinap = 0 and limp’cosap=0, p> O0,a 49. 
px ps 

For, if d(a) is any test-function, 

1 


A x) p? sin xp da p? 


-4(x)cos xp : Rew oa 1 d'(x)cos We) dx 


p? 1 d’(x)cos vp da 


L 


and 


| A( x)p?COS vp da p* 1\ d/( y)sin xp ; ’ p” 1 ( d’ ( v)sin vp dx 


£ 


p” I dh’ ( v)sin xp da. 


A repeated application of this process leads to the inequalities 
d(x)p” sin xp dx| : p l AP (x)| da - Kp 1 


x f 


and | d(x)p” cosap da} < p™ ( AY *V(x)) dx < Kp-', 


x 





x 





so that 
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» that the limits as p tends to infinity are zero. Consequently, in the sense 
fthe theory of distributions, by repeated integration by parts 
R 
: ee 

: ). s1n ap 
lim p” sin xp dp lim ( \)p! aa 

I yl? J CcOosSap 

0 


1U, if x = @, 


Ip 


\oo0 ~=6owwith pole of order p+-1, if « 0. 


This suggests that 
lim | p”’sinapdp = K &(a), 
j . 


ere 8”)(x) is the pth derivative of the Dirac measure; and 


9 , ‘“ r du 
lim H(t)l(R) kK (1 —u?)2?4-BD54—2)(ct—ru) 


choosing A to be real 


constant AK is determined by the inverse relation 


K (ce \*-* (ct? , sin 2zcpt 
r=) | 1 | c f ct,t > 0), 
[ - 27 p 

nd after a somewhat tedious calculation involving Lommel’s function 

11, 374, 345) it is found to be 

| 
wc" 3( fo 2)! 

here w, is the area of the unit hypersphere. The elementary solution for 

31s thus 

| 1/e\*-2/e2#? (k—3 
5) (G 1] H(ct—r). 
D.C" (A 2 ct] er 

). The initial value problems 

Che solutions of the initial value problems may now be found by means 
{the convolution theorem (iv) of (4). In the initial value problems 


B u,(x) d'(t) u(x) d(t) | c?, 
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as may be seen from the Gaussian integral theorem for any self-adjoint 
differential operator of the second order (ef. 3, 435) 


| en } (dL 


fy « 


(V) 





w| u L\¢ |) dx, day...dx;,dz 


¥ ¥ i ou e c db OX; 
| nas 6» a,.— —U S Asp. ‘da 
J’ —_ — : 


p. Cc Xi, c Xp, CV 
(S) 


in which the given values u,(x) of uw are associated with the transversal 
derivatives of don S and the given values wu, (2) of the transversal derivatives 
du/dt are associated with ¢ itself. 


When k = 1, 


E (x,t) = —.H(ct—|z\), 
x r+c(t—r 
# x)d'(t a ae F ae 
E(x, t) * ol) 0 ( ) —= 0'(r) dr | uy (E) dé 
Gc aC 
x r—c(t—T) 


; a 
and Ele, t)* te : L oe u,(&) dé 
c? a a a 
r ~ct 
u,() dé, 
2c J 


so that we recover d’Alembert’s solution (1). 


When k = 2, c H(ct—r) 


E(x, t) “ae ; (cP 72) : 
E(x, t) * w(x) 3'(t) s vias H{e(t =*} p]uol2—6) dédy 
c re J 4 Ner(t—7)*— p*; 
— — —- { 8'(7) Uq--)(U9) dr = -— © pte) 
where p? = (x—£&)?+(y—7)?, and 
E(x, t) x a Me(Uy); 


c* 27 


so that we recover the Poisson—Parseval solution (3), section 1. 
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When A 3 


E(x. t ( d(ct Yr) 
Aor 7 
: l)(.v) 0 (7) l F P io € i dl c(t T) 2 9 
E(x.t) x , (7) dr | | : f luo()p2 dpdQ, 
c* ATC Ida p 
ert 2ti(y (z—C)? and dQ is the element of solid angle, 


| 0'(r) dr c(t T)u,(€) dQ 


wo « 
r=c|t—T 


t 7\(Uo) dr 


( 


6’ (7) 4are(t—7)M 


1 
d 
- V ly) 
dT 
U,(x) o(t 
nd similarly E - at tM,,(u,). 
72 
When 3 
| | k-—2 ct? i(k —3 
I (=) | | H(ct—r) 
W 2)! / ot i 
I r)5'(t l¢ k-—3 . _ Vu j 
ae (c2t2— p2)!&-B pM (up) ¢ 
; aia) | p p ol o) P 
0 
in easy ileulation l 
ct 
%,(2%) O(T) | o\*-2 ¢ 49 9 ‘ 
E(x, t) x i | (ct? — p?)'*-)p M(u,) dp, 
c* ( A 2)'\ot 
0 
that 
Z {P cf? r) (hk Xr Me (ug) dr 
k—2)!l\o8 
OP 0 ce se | 
| | (c?#2?—1r?)& Dp M(u,) dr}, 
ot J 
0 
veneralization of (4), section 1. 
6. The inhomogeneous problem 
The solution of the inhomogeneous problem is given by the formula 


u(.v, t) Exg 
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which leads immediately to the formulae (5), (6), and (7), and the following 


generalization of (8) 


( 


\kK-2 7 
| (c2t2—r2)§4-Dy Me fg(r, t)} dr. 


« 
0 
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PROPAGATION IN A FERRITE-FILLED WAVEGUIDE 
By LL. G. CHAMBERS (University College of North Wales) 


ved 2 December 1954] 


SUMMARY 


Perturbation methods a1 sed for the solution of the fields in a waveguide filled 

ferrite material which is subjected to a static magnetic field in the direction 
its axis. It is shown that quasi TE and quasi TM modes exist and the first terms 
the expansion are calculated for the case of a waveguide of rectangular cross- 


1. Introduction 


Ir has been found experimentally that, when a static magnetic field is 
impressed on materials such as ferrites, they become anisotropic. The 
nechanism by which this happens is not altogether known but a possible 
explanation has been given by Polder (1). The atoms comprising such 
substances are arranged in a cubic lattice and behave as magnetic tops, 
heir angular momentum and magnetic moment being parallel. If an 
external field is applied the tops will align themselves in certain directions, 
nd if small disturbances such as a microwave field are superimposed upon 
he static magnetic field, the additional magnetic induction and the addi- 
ional magnetic field intensity due to this microwave field obey the following 
relation 


B nH Je HL, Sn l 
Bb, ju H, T pH, 
, = *,.o. (1.1) 


\ time factor expfiwt! is assumed and the static magnetic field is in the 
- direction. The exact manner in which yp, yp’, », behave with this applied 
static magnetic field is not important from the point of view of this paper, 
ind has been given by Polder (1) and Van Trier (3). It may be noted, 
however, that («’/u) is proportional to this applied field. 

Ifa waveguide filled with ferrite be subjected to a static magnetic field 
n the direction of its axis, it is possible by using Maxwell’s equations to 
determine the transverse components of the electric and magnetic inten- 
sities in terms of the longitudinal components. Without going into details 
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of the analysis which is given elsewhere (2), the field is given by the 
equations 
(x*—k"*)E, = V(xk*yF,+-wy*' H,)- 
9 i, x Viw(K2u— kh’? py’) H,— yk’? EB, (1.2 
(«3 k’*)H, V(K yH., wk’ E.) 

Ji, ™ Vikrwek,+-yk'H.}, (1.3 
where a factor exp(yz) is assumed, i, is a unit vector along the z direction 
E.,, H, are the field components in the z direction, E,, H, are the field com- 
ponents in the transverse directions, and 

K wep -+-y?, ,"* wen’, (1.4 
where e€ is the dielectric constant, V is the planar gradient vector. The 
values of FE. and H. are determined by means of the two differential 


equations 


V7E, + fn? —(k'2y’)/p! E,—wyp(p' |p), (, (1.5a 
V?H_--K?(u,/n) HH 4+-wye(p’ py) BE, = 0, (1.5b 
subject to the boundary conditions 
* E 0, (1.6a 
, jay ee cay oe jolt — btn’) = 0, (1.6b 
on cs on 


where €/cn and @/és represent the normal and tangential derivatives on 
the boundary. 

Because of the nature of the boundary condition (1.6b) it has not yet 
proved possible to solve these equations generally. Solutions only appeat 
to be possible for waveguide cross-sections, such as the circle, for which it 
is possible to express an are distance along the boundary in terms of one 
coordinate only. 

However, if the ratio (u’/uw) is small, it is possible to use perturbation 
methods, and a solution may be found in terms of a power series in (j' 1 


The first terms of the series will now be found. 


2. Preliminary analysis 

It will be observed that in equations (1.5) and (1.6) the following 
quantities may be expanded in terms of (u’/u), namely H,, E,, y. Now) 
occurs both as y and as y*, and it will accordingly be convenient to recast 


the equations so that y?, rather than y, is the unknown quantity. There: 


Equations involving the boundary conditions will all be denoted in this manner \ 
other equations are applicable throughout the medium. 
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fore, writing wuH QV, yE Pr, wwe, a p/p, and ¢ He/ fs 
equations (1.5) become 


V?P+hy?+kh?(1—a?) P+lay*Q 0, (2.la) 
V20-+{y?+k3iQ—kaP = 0, (2.1 b) 
subject to the boundary conditions on C, 

* P = Q, (2.2 a) 
‘ X 208 ) = k?(1—a?) 1 OW 0 (2.2b) 

Os On on 
[t will be seen that, in this form, y occurs only in the form y?. Now if « = 0, 

equations (1.2), (1.3) become 

1 E, = VyE,+ji, x V(wpl,), (2.3) 
«*H VyH,—ji, x V(weFk,), (2.4) 


the cross terra in equation (2.1) disappears, and the fields split up into those 
for which # 0 (TE) and those for which H. 0 (TM). Moreover the 
modes which occur in the ferrite-filled guides may be expected to be of two 
series, one of which reduces to the TE type, and the other to the TM type 
a when « = 0. Further, to a good degree of approximation (3), € = (u,/y) = 1 
und, as the interest lies jn the behaviour with a, this value for ¢ will be 
ssumed (the extension to € + | will be almost obvious), whence the equa- 


Tions hex ome 


V2P 2+-k?(1—a*)!P+ay?Q = 0, (2.5a) 
V2 +(y?+kh?)Q—ak?P = 0. (2.5b) 
D { . . . - 
Before proceeding to the perturbation, a further mathematical tool is 
r juired—the Green’s functions G, G® which are defined as follows: 
V2qQi.2 k?+-5 2) {0,2 d(a— 2" )d(y y’), (2.6 a) 
Ge 
Cur 
G Q), O ond, (2.6b) 
Ch 


When these are known, the solution of 


V?y+-(k?+-yo)x = f(x,y), (2.7) 
ibject to the appropriate boundary condition, is given by 
iV ie . 
v,y | f(a’, y G(x, a's y,y') da'dy’, (2.8) 
7 Where the integration is over the cross-section of the waveguide. The 
Green’s functions G™, G® are calculated in Appendix A for a waveguide 


ot rectangular cross-section. 


































438 





Li. G. © 


3. Calculation of perturbed field 


the form 





HAMBERS 


The field, and propagation constant, are assumed to be expressible j 


P Y a"P ». 
n=0 
0) Y a''Q, 3.2 
n=0 
2 > ytys, 5.3 
n ) 
where y? is not necessarily positive. We now substitute for P. Q and y2 
in equations (2.5) and (2.2) and obtain the following formulae: 
Equation (2.5a) becomes 
> "(VAP +kP )—a®k? S&S a" P. 
? 0 n=0 
2.@¥%. 2 OF +e SF ays & a"Q, = 0 
n 0 ; 0 _=0 i 
or 
n 
Y atl (V2+k2)P + +f... 
n=0 l r=0 
} n 1 
2 > ee Ya" ¥ ¥Q, : 0 
n=2 , Pres r—0 
whence 
V"F, (hk? va) Py 0, (3.40 
V°P, + (k?+-y2)P, v2 P,— v2 Qo, (3.41 
n 1 n 1 
2p a 2 . 2 2 ‘ 
V2P +-(k Yo)P oe gl > ¥P._.—EP,_. 3 G1 
r=] P=0 
. (3.42 
Equation (2.5 b) becomes 
> x" (V2+-k2)Q,, > yay? > xv", —ak? ¥ a" P 0 
n 0 n 0 ” 0 n 0 
or > a«"|(V2+k2)Q,, > v2Q9,--|—k2 ¥ a" P 0, 
n=0 r= § an 
whence ‘ . 9 
V°Q (k? V0) Vo = 0, (3.00 
V°Q, (k? T v0); vi Vo k*P,, (+ 5] 
. ; ne n J 
V2Q,, + (k?+y2)Q, M— > v¥?Q, k2P _, (3.52 
7 ‘ f=] 
Equation (2.2a) becomes 
m Y a'P 7) 
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wnence 


P 0 (3.62) 
juation (2 2? b) becomes 
CO oP 
* % \ : jk*a S x! , 
Pan hoi Cs ead On 
0 n 0 


CO 4 7 : a) 
ea S 2S nna gS an ¥ pn 
a — Ss — cn ‘ — — ar con 
? l n 0 T 0 
co . * ra) 
jee S a jee S om te? = 0 
— on — ch 
0 n 2 
 eQ 7 
i(k? y2) v0 Q, (3.70) 
Chl 
() P », OW . «oY 
2 t j(I yz) v1 jy; vo 0 ( } 71) 
( } M ch 
() / ; ,oV aa) .. ! 
* i] j * 2 Cn —7 jk? v) pec en-2 0 
homed — cn Cn ch 
0 
(3.72) 
It will be observed that equations (3.42), (3.52), (3.62), (3.72) determine 
Q) n terms of the boundary conditions and the lower P,, Q,,, v3. 
Thus from P,, Q, it will be possible to build up the quantities P,, Q,. 
The Té ot zero orde are 
V?2P,+(k?+y2)P, = 0, (3.40) 
v 7), (}:2 V5) Qo 0. (3.50) 
P, = 0, (3.60) 
2, OY 97 
j(k?+5) —" = 0, (3.70) 
Ch 
, ) 
which 1s eg ivalent to 4) 0, 
Chl 


These represent the TE and TM modes for a guide filled with a material 
whose electrical constants are e, 1 but with p’ 0. (These will, of course, 
be the same in their cross-sectional behaviour as those for an empty 
guide.) On examination of the equations (3.4)—-(3.7) it will be seen that 


ny he 


ferrite filled guide will be composed of those modes such that 
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P, = 0 and those modes such that Q) = 0. Thus the solution of equations 
(3.4)-(3.7) may be considered separately for P, = 0 and for Y, = 0. Modes 
such that P, = 0 will be termed quasi TE and those for which Q, = 0, 
termed quasi TM. The first type will be dealt with in section 4 and the 
second type in section 5. (Equations in these sections which have appeared 
previously have their previous number on the left.) 


4. Quasi TE modes (P, = 0) 
The governing equations in this case are 

(3.50) [V2+ (k?-+-y2)]1Q, = 9, 

OM 


on 


(4.] a) 


*(3.70) 0. (4.1b) 
Q, will be taken to be any particular mode of the empty guide given by 
equation (4.1). 

The first-order terms become 


(3.41) V2P + (k?+-y6)P, ye Q); (4.2a) 

(3.51) VQ, + (k?+y72)Q, V3 Qo, (4.2b) 

*(3.61) F, = 0, (4.2¢ 
" ee 

mia th) v2 Xo _ y(k?- y2) © : 0. (4.2d 
cs on 


From equations (4.2a) and (4.2¢) it is possible, using the Green’s function 
G® to calculate P,. This is done for a rectangular guide in Appendix B 1. 
There remain to be calculated yj and Q,. For the calculation of y? the 


method of Caprioli may be used (4). From equations (4.1 a) and (4.2b) 


vi | Qod8., 


| (Q,V?Q, 


© 0 


0; V?Qo) dS 
S 


and by Green’s theorem 


, Q Q f Q 
(Q, V7Q, 0) V2Qp) dS | (20: =I 0; = "| ds @ “ad ds 
Cc 


S C4 


S being the cross-section of the guide and C of the wall. This follows 


) ¢ 2 : - 
. ($Q5) and we integrate completely round C. Thus 
os ° 


oq 
because Q,— 


cs 


v2 | Qe dS QO. 
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ind so yi = 0. (4.3) 
Equation (4 2b) becomes 
V2Q,+(k?+y78)Q, = 0. (4.4) 


We have by Green s theorem 


tl 


| [G42V2Q,—Q, V2G®] ds’ 


on n 





Y Go 
Gre _g eG | as (4.5) 


Cc 


where S is the cross-section of the waveguide and C is the wall. Using 


equations (4.4), (2.6a), (2.6b) it follows that 


V(x ,y jo x )oly—y ) dx dy Gi2)_¥1 ds’, 


hencs 
cYV 
) Yu Ge t1 ds’ 
Ch 
v= 4 ' ; cy (a y’) : 
Ge. 2": 9. 7) (4.6) 
- . cs 

his procedure is carried out for a rectangular guide in Appendix B 2. 

Thus P,, Q,, yz have all been determined and it is possible to go on to the 
rms of | onel order. 

Now equations (3.42 3.52), (3.62), (3.72) may be regarded as ex- 
ressing the P,, Y,,, y= in terms of the quantities of lower order. Remem- 
bering that P 0, equation (3.42) may be expressed as 

) 2)] 7 " 
\ k?+-y5)|P,, P.2 (4.7) 
je 1 
*(3 62 P 0. 
r,_, being terms of order not greater than n—1, and the solution follows 


mmediately by using the first Green’s function as was done for P,. The 


valuation of und Y,, presents a little more difficulty. 
Using Capriolis method again. 
r cY , ow) oY . i ° 9 ’ 
Qo en ds — | (Q,%2"_@ <9) ds , | (Q,V2Q,,—Q,, V2Q,) dS, 
) Co? Ch 


S 
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which may be rewritten, using (3.52) and (3.72) as 


n-1 





= | Q, \ 2 0@, 4 jae 1 
i(k? + Yo) A oa Cs en 
1 : 
> ( “7.9 ( 9 
j * yp nar 4 je @n-a] dg 
hat on on 
r=0 
| I? Os T k*P, 1 Vor Q > “, Yr 0, | ds, (4.8 


whence y? may be calculated. 


To obtain Q,, write equation (3.72) in the form 


< 20Q a xX ,0Q eQ eQ 
: 2 n-1-r 2 n-1 2 n—-r 1,2 n—-2 */ 7,9 ’ 
+ > #- Lh i> 4+ - j(k?4 Yo) 
— cs on = on on n 
n 0 r 


(4.9 


The left-hand side is known completely and so it follows in exactly the 


same way as equation (4.6) was obtained that 


Q,,(2.y) G2) ds’. (4.10 


5. Quasi TM modes (Q, = 0) 


The fields in the quasi TM modes may be calculated in exactly the same 


way as those for quasi TE modes. It is, however, of interest to consider the 


first-order terms 


(3.41) V*F | ( f:2 vs) P; yz Fos (5.1la 

*(3.61) P, 0. (5.1b 

(3.51) V0, 1 (fy? v5)Q, a (5,1¢ 
‘ 7 7, o¢ ‘ 

*(3.71) peo j(k? v2): 1 0. (5.1d 
: on ' On 


Applying Caprioli’s method, 


0= | (ASRS ds — [ (AVR-RVR) As 
m Cn Cn 


whence 


It follows from equations (3.41) and (2.8) that P, 0. 
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n exactly the same way as was indicated previously in section 4 
” rye) 9@s 3.0 
O(a, y | G2)_*1 ds 
J Con 
. : ; 
YO 2) OL 0 407 a 
YO . | G2 0 ds’. (5.3) 
k? YO. on 
“ 
ilculation of this quantity for a waveguide of rectangular cross-section 
s given in Appendix B 3 
The higher order terms (i.e. 7 2) may be calculated in the manner in- 
cated in section 4. 


6. Discussion 


The transverse fields may be calculated from equations (1.2) and (1.3), 
) the field is completely known. It will be observed that the terms of 
rder a in the expansion of y? vanish and thus, because a = (u’/) and 

is proportional to the applied magnetic field (3), the propagation con- 
stants in a waveguide filled with ferrites are insensitive to the effect of a 
iall ma etic field 

The utility of this lies in the fact that once the quantities P,, Q,,, yz have 
en calculated, it is possible to find the fields for any arbitrary value of «. 
The quantities P, and Q, for the quasi TE case (P, = 0) and the quantity Q, 
for the quasi TM case (Q, ()) have been calculated in the appendices for 
e case of a waveguide of rectangular cross-section. For this, an exact 
lution is difficult to obtain as a boundary condition of the form 

seg pe 
Cs Ch cn 

uses trouble, the difficulty arising from the fact that over part of the 
boundary x is in the x direction and over part in the y direction, and 
similarly for s. This difficulty does not arise in the case of a circular 


oundary where s is always in the @ direction and n is always in the r 


urection 
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APPENDIX A 


Calculation of Green’s functions for a rectangular guide 


Let us suppose that the guide is given by 0 x<a,0 y b. 
ple gland 
: : 71 7 
Let Ki, = yy (Al 
‘ sf a> 62 _— 
€ l m 0, 
1 49 
2 m 0 re 


The Green’s function of the second kind is given by 


eG) eGge oe rn ao : 
= ~- ———— +} on GI o(xr—azx )o(y—y ) (A3 
Ci” Cc y 
and 
oG 2 
. O over a 0, x a:y O, y b. (A4 
on , : 


It follows from the boundary condition that G@) must be a double Fourier series of 
cosines. Also 
G(r,’ sy.y’) GO)(a’,x3y',y). (A5 


Accordingly we write 


Ge 2 ’ Ys y’) : z Zz Gp €, cos(= Jeos(“™ -Jeos()cos(=4 : A6 


0 s=0 


which satisfies all these conditions. Also 





s r 
Yale) 22) # ‘ , 
c2Ge 7G = ¥ oi VTL VTL S7TY STTY 
—-- -—" TK Gi S S €, €. G - K2, cos| - Jeos| -)eos| “)eos| }, 
ox* cy” y ap : a a »/ ) 
0 s=0 
Al 
; hoe : ‘ l mT ‘niy’\ . 
Multiplying both sides of the equation (A 3) by cos| )eos| _ -), integrating 
: * ab a b 
over 0 2" a, 0 y’ b, and using the result 
a 
Lf "Ta mare’ 
€,COS COs da oO AS 
a a a 
0 


it follows that 


x 


ab S S G K ia 5., coal Jeos(*™) cos(—— leos(“7), \9 
pa ° a b a 5 

T Us 0 
whence G (abx?.)7} Al 
and 

G(r 2a’ yy’) >. * e,€,(abn2,) conf bonah o bounl te benef ). 
a kod hee a a b b 
os 0 


In exactly the same way it may be proved that 


x x 
GM (a, a soe") S S 4$(abn>,) 1 sin(")sin(™ ) sin (*74) sin(=™ ). 
a =) a a b b 


ls=1 
(A 12 


This clearly vanishes on the boundary. 
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jalculation of }, for quasi TE mode (Fe 0) 


— a 
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Koo Fy Ys Qe (B 1.1) 
on the boundary. (B 1.2) 
P 1¥ Qo(x’, y’)} da’dy’. (B 1.3) 
rt} n by 
WATT A MTTY 
) ys Cos (B 1.4) 
(——]eos(—=) 
} N NS { } 
a ~ —T STTY stry’ mora’ nory’ . . 
n{ } }sin{ sin “~Jcos Jcos(- *) da dy 
| | b | 5 a b } 
> % P / n| }sin( } (B 1.5) 
2 | co mar sin ral dC. (B 1.6) 
alculation of (, for quasi TE mode (P, 0) 
CUo(r .y) ’ 
O,(2x, |G r,x’sy,y’) ~~ ds’. (B 2.1) 
} cs 
r tl tion, 
r 0 ,©Y(a.y’) 
1 ( ( da G@)(2, a; y, y’) ——_ dy 
‘ . cy 
0 
0 
Osx" 6). . 7 cQ, (0, y’) 
° a as G'2)(x, 0; y, y’) ———— dy’, 
Ca f itn cy 
4 
r integ1 1 iting integration in the anticlockwise sense along 
wall ry ve rewritten in the form 


(B 2.2) 
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where 
a 
- TT oY, (x , 9) ( Qo(a’,b) 
EB e COs } = COS S87T - da 
a cu CA 
0 
bh 
stry’ CVO, y’) oQ, (a, y’) 
cost } —- COS?TT > da 
. P cy cu 
0 
Again using the fact that a typical value for Q, is 
. co 
mary niry 
cos| : Jeos| = ), 
a b 
a 
: rare’\. mar}. (mare ' 
E.. | cos| }| | -coss7 cos nt | sin ) da 
a a a 
0 
b 
sTry’ nT nqry 
| cos| }| l COS 177 COS TT | sin| ) dy 
b b h 
0 
qT nT, 
> [1 +-coss7cos nar jl,» ~ {1 cos rm cos ma |l.,,, 
whence 
¥ . 
Ki, Q hy ‘ * Vara STTY 
; = S S €, €,(abK?.) Leos Jeos| : 
Vo yay a a b 
0s=0 
(mr nq 
‘ie [1 COS 87T COS nat |e r > [1 COS ITT COS INT I \* 


3. Calculation of (, for quasi TM mode (Q, = 0) 


ye . . a: Pia .yu) . ne 
Q(z, y) i Galax Pe | ) = _— ds . (B3. 
Koo. on 
whence 
a 
Koo —_ CP(2’,0), , 
> (x,y) | G)(x, 2's; y, 0) stale dx 
Yo M (al 
0 
: ,, oP, (a, y’) 
| G?)(a,asy,y’) ~ y dy 
0 
0 0 
q eP,(a’.t i Crete.) 5, 
G'?)(a, 2’; y, b) ia de — | G(x, Os y,y ) of z dy 
4 cy 2 ° Cx 
a b 


each of the four integrals representing integration in the anti-clockwise sense along 
one of the four walls. 


This may be rewritten in the following form: 


a Q1(%,Y) > Zz €,€,.(abK?.) 'cos(")cos("F") Fre (B3.2 





where 


Again 


when¢ 
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F stry’\ [oP,(0, y’) oP,(a,y’) , 
| cos(—" ) COS rt dy’. 
b Ca 


CA 


al value for P, is 
pom NITY 
sin } sin( : ), 
a b 


COS 87T COS 117T da 


’ stry’ niry \mt x 
cos| )sin{ ) | 1 —cosra cos mr] dy’, 
b b a , 
wthgT b 
77 } =) cos ri cos m7}, 
= a 7 


mar b | 
COS ITT COS mi }|. 
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TWO PROPERTIES OF SPHERICAL HARMONICS 


By HAROLD JEFFREYS (St. John’s Colle ge. Cambridge ) 
[Received 27 January 1955} 


SUMMARY 


The integral of the square of the gradient of a solid harmonic over a sphet 
evaluated ; the corresponding integral for the second derivatives is also evaluat 
and the results are applied to an integral that includes the elastic energy in a 
strained sphere and the rate of dissipation in a viscous sphere. 

A natural definition of the irregularity of a function over a sphere leads to t] 
conclusion that the irregularity is stationary for small variations of the functio 
when the function is a surface harmonic and that the irregularity of any functior 


is greater than that of the lowest term in its expansion in surface harmonics. 


THE first of the properties in question has been used in geophysics (1, 2 
but may have other applications. The first application given here is used 
in a later paper. 
1. Let A, be a solid harmonic of positive integral degree n, and put 
kK, rs... (1 


We show that over a sphere r = a 


ea: ae n(2n- 2 46 
| ‘| dS ( an) A: ds, (2 
Cx. 
We have by Green’s theorem 
ok, \? 2, OK ; n Siena - 
| | (- | dr = | ‘x "46 | | K* ds (3 
CX; r CX; a 
. e e i « e i e e 
Per yn (aK \2 
) ok 
| | 7 drdw 
2n-—2 * 
a a C2 i/r=a 


2n- = | | é& ‘) ea " 


whence the result follows, since dS a*dw. 


and also 


1 


double summation (1, p. 74) 


{ { he a =) as (n—1)n(2n—1)(2n+ 1) [ [ K? dS. (5 


ry 
; OX, a 


Since @K,,/éx; is a solid harmonic of degree n—1, repetition gives tht 


(Quart. Journ. Mech, and Applied Math., Vol. VIII, Pt. 4 (1955)] 
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1.1. Let uw; be a tensor of the form (1, p. 81) 


o*k ' Pala oK,, 2," K,, 


CX, OX CX, Ca 


+- Fi, 2; 2, .K, + F,8;,K,,; (6) 


U 
t/ 


where F, F; F.. F, are Tunctions of ,. Then 


(nF,+rF,+3K)K,,. 


— 
~I 
~— 


The sealar 
oe r2| e2k ) "(r, ' 2(n -1)F,\? 2(n = ries a 
: : 12 v, OX, 12 i r a | \ Ox; 
{1 n2 F's ‘nr, F; t+-2n(n 1)F, F, T gr? F'2\ K2 
(/. eK, \2 (.aK,\2 _ 
Geli r2 2(n—1)*|2 Ln?(n—1)#K2} + 
(? Ox. C2 | \ Sy \ ) 
' kK asm oe 
103|( - "| n2K?| +2G5K", (9) 
= Ox, | 
where 


Then the integral over a sphere of radius r is 
| UW: us». aS iIn(n l)\in+1)(n- 2)G?4 n(n 1)G3+ : 


(11) 


If u,.is a stress, and thus w’,. w;,. the Mises function, the G, terms represent 
the contribution from the hoop stresses, G, from the shear stress over 
horizontal planes, and G, from the difference between the radial stress and 
the mean of the greatest and least hoop stresses. 


If u-is a displacement, and 


the strain component ¢ S 
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and we find, with w,, Cun 
F F’ I 
Gy ; G, 1 Gt 2(n L)— 
"Sis al if an 
, i onin—1).. 
“3 - n| a] F,. (14 
= r y+ r 


This way of expressing the strain components is much more convenient 
than those that introduce derivatives of: K,,/r?"*1, It seems to have bee 
introduced by Love (3, p. 15). Another type of solution arises when there 


is tangential stress over the outer boundary (2). 


2. Inthe same conditions, if 6, y are polar coordinates, and f is a functior 
of 6, yx, a natural estimate summarizing the general magnitude of f is 
|| f*dw. A natural estimate summarizing the rate of variation of [ wit! 


position on the sphere is 
"re (lof\? (1 af\ 
I i {|| (— Py | dev 
JJ \e6) sin cy} | 
Then the ratio of these two expressions is a precisely defined quantity that 
can be regarded as a measure of the irregularity of f. It is unaltered by 
change of axes. We have for any function of position 


(ed\2 l 2éd\2 ( (ed \2 d\2 
ia ie a a (| (" | , | 


/ 
and if ¥ Ss, d “ pn nh 2 


9 


( =a | dN 


a 


n(n 


a 


WT on 7a : 
: | | Kk? dw n(n 1) | | fP dw. 3 


replaced by f+-f. 
i {ff (lef\? 1 ef\? a. 
si | | : a} a A =| —Af . sin 6 dédy 


In 
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. on 


dé 
x 0 


| ns 
+ ss 
sin @ ox 





; Gaal” 
Ss Gg z of ly 
inde, ax | 
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ee 
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f l of 


Cc ‘ c 
(sin é ; : 
o@ fal] sin 6 cx 





Asin Of |} dy 





Further, by the methods of the calculus of variations, we find that if fis 
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[he integrated parts vanish and the condition that the integral is stationary 
for genel il variations of is 

ef l ef 

ising ee ~~. 

sin 6 c@ rat sin“@ oy" 


Af Q). (5) 


which is the differential equation satisfied by S,, provided that A n(n-+1). 
Then. bythe usual argument, our estimate of theirregularity of fisstationary 
for small variations of f provided f is a surface harmonic of integral degree n, 


] 


nd the estimate of irregularity is then n(n-+- 1). 


It is sometimes said that spherical harmonies are an artificial method 


inalysis, and therefore it seems desirable to point out that they have 
direct relation to the property of smoothness when this is defined in what 
seems to be the most natural way. 


It n be added that if A,,, A, are any orthogonal harmonics, 


ft Ok, Om. 


dS — 0. (6) 

Cu; CX 
us fol general f there are theorems similar to those that hold for 
spansions in orthogonal functions. A consequence is that if the expansion 


fin harmonies contains no harmonics of degrees less than n, the irregu- 
ty of fis at least n(n 1 

\ problem in geodesy is the interpolation of gravity from irregularly 
stributed observed values. This result shows that the smoothest inter- 
late is that given by fitting surface harmonics up to a given degree; this 


so takes account of the fact that five of the lower harmonics are theoreti 
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NOTE ON THE IMPROVEMENT OF APPROXIMATE 
LATENT ROOTS AND MODAL COLUMNS OF A 
SYMMETRICAL MATRIX 


By W. J. GOODEY (62, Welland Lodge Rd., Cheltenham, Glos.) 
[Received 21 October 1954] 


SUMMARY 
The method given below for the improvement of an approximate set of latent 
roots and modal columns of a symmetrical matrix is essentially the same as Jahn’s 
second method described in (1), but the proof is here set out entirely in matrix 
notation. 
The present paper is, in fact, the counterpart of the notes by Collar (2) relating 
to Jahn’s first method. 


1. Ir VM is a modal matrix of a symmetrical matrix A, both M’AWM and 
M’M are diagonals. Thus, if m is an approximation to M, m’Am and mn 
will be approximately diagonal, and we can make use of these products to 
obtain a better approximation to M, as follows, 

Let m’Am = A,+-aandm’m = A,+5, where A, and A, are the diagonals, 
and a and 6 the non-diagonal parts, of these matrices. Further, let 
M = m(I-+ pn). By suitable post-multiplication by a diagonal matrix we 
can always ensure that /+-y has unit diagonal elements. Since such post 
multiplication of W leaves it still a modal matrix of A, we can suppos 
without loss of generality that » has no diagonal elements. 

Then (J+ p’)(A,+a)(1+p) and (I+-p’)(A,+5)(1+p) are diagonals. 

Expanding these products, and neglecting second- and third-order terms 
we see that A,+-(u’A,+A,u+a) and A,+(p’A,+A, u+56) must be diagona 
matrices. But the items in the brackets have no diagonal terms, and s 


WA, +A ,pta = p’A,+A,p+b = 0, 
whence A, pA,—A,pA, = bA,—aAg. (1 


If the term in the rth row and the sth column of p is p,,, the corresponding 
term in (A,uA,—A,pA,) is (AyAo.—Ao-Ars)H 
are respectively the elements of A, and A,. Hence p,, is derived from thi 


re Where the A,’s and A,s 
corresponding term in (bA,—aA,) by dividing by (A,,,A,,,—A2,-A;,5): 
Since the latent roots of A are A,/A, (with errors of order y*), the method 
will not work if two or more roots are equal, since the corresponding 
divisors then vanish, and it is clearly not very effective if any of the latent 
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roots are close together. In general, however, the method is a very attrac- 


tive one, inasmuch as it produces an all-round improvement of the modal 
natrix in one operation. 
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VIBRATION OF CERTAIN SQUARE PLATES HAVING 
SIMILAR ADJACENT EDGES 


By HUGH L. COX 
(Research Engineer, T. ¥. Lin and Associates, Los Angeles, California) 
Received 17 February 1955) 
SUMMARY 

The fundamental frequencies of flexural vibration are determined for thin uniforr 
elastic square plates that have two adjacent edges either clamped or simply supported 
and that always have the other two adjacent edges free; a pin-point support exists 
at the intersection of the free edges. Finite difference methods, which simplify tl 
treatment of the free boundaries for definite values of Poisson’s ratio, are used iy 


conjunction with extrapolation procedures to obtain the approximate solutions. 


AN approximate solution is presented for the fundamental frequency of 
flexural vibration of a uniform elastic square plate that has two adjacent 
edges clamped and the other two adjacent edges free; the plate is point 
supported at the intersection of the free edges so that this point cannot 
deflect. The fundamental frequency of an otherwise similar square plate 
which has simply supported adjacent edges is obtained also. (See Fig. 1. 
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VIBRATION OF CERTAIN SQUARE PLATES 455 


[he solutions are valid for plates which have a value of Poisson's ratio 


£ (2 
The usual assumptions regarding the bending of thin elastic plates are 
ude. If the static lateral loads on the plate are replaced with equivalent 

nertia forces when the plate is vibrating freely, the governing differential 

equation Tol the plate is 
= —_ qu 
V2V2w0 ! = 0, (1) 
gD 
where s the plate deflexion, ¢ is the lateral weight per length squared 


{the plate, D is the plate stiffness, g is gravitational acceleration, and 
is natural frequency in radians per second. Equation (1) may be 
written in finite difference form at a selected number of points. First-order 
perators, which include Poisson’s ratio, for points on and adjacent to the 
free boundaries have been tabulated by Jensen (1). If equation (1) is 
written in difference form, the resulting simultaneous linear algebraic 
equations may be written in matrix form as 
f*A—rAI|W 0, (2) 
where f is the plate width 6 divided by the grid spacing, A is a matrix of 
coefficients of the V?V2w terms in (1), J is a unit matrix, W is a column 
natrix of deflexions, and 
qui b4 
gD 
Equation (2) is in standard eigenvalue form and the lowest root A may be 


letermined rapidly on an electronic digital computer. Equation (3) may 
pial | . 


ne rewritten as 
k (gD\3 
Wy = 5 (4) 
h2 q 
where \'. For the case of a square plate having two adjacent edges 
clamped, values of / have been found to be 
i 12-89906 for f = 6, 


12°55425 for f 5. 


Values of / for other values of f indicate that k is approaching the true 


ilue monotonically; consequently, extrapolation (2) vields 
6- 5? oe ” - 
} 12-89906 12-55425 13-683. (5) 
62— 5? 62— 5° 

= 13-683 /qD\4 

Thus w ; (" )*. (6) 
l h2 
a: q 


The above result probably deviates from the true answer by less than one 
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per cent. Extrapolation gives very accurate results (2). The fundamental 
mode shape is symmetric about a diagonal line passing through the point 
supported corner. If the pin support at x = y= 6 is removed, the 
coefficient k for (4) has been found by Young (3) to be 6-958. It is seen 
that the post support at the corner x = y = b roughly doubles the funda- 
mental frequency of a square plate having two adjacent edges clamped and 
the other edges free. It is interesting to note that the second mode of the 
plate in Fig. 1 will have the same frequency as the plate discussed by 
Young since both plates have a nodal line along the line 2 = y for the 
second mode. The value of k for the second mode (3) is 24-08. 
For the case of a square plate having two adjacent simply supported 

edges, values of k have been found to be 

k = 8501848 for f= 5, 

k = 8222674 for f = 4. 
Extrapolation yields 
ne 42 
k= — , 8-501848 — ——___ 8-222674 8-998. (7) 
ra 5s 4° 
se 


Thus, wy, 
b2 q ; 
\ / 


(8 
The above answer probably deviates from the true answer by less than 
one per cent. 

If the plate in Fig. 1 rests on a uniform elastic foundation (4), 


13-683 /gD\4  (gK\3 
Wy = : “+1 {° . (9) 
b? q \ q 
for the plate having two adjacent clamped edges and 
8-998 /gD\}  (gK\4 
Wy - = : “+t f° Pe (10) 
O° \q q 


for the plate having two adjacent simply supported edges. The symbol A 
represents the elastic foundation modulus per length cubed. 
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TABLES OF BENNETT FUNCTIONS FOR THE 
TWO-FREQUENCY MODULATION PRODUCT PROBLEM 
FOR THE HALF-WAVE SQUARE-LAW RECTIFIER 


By R. L. STERNBERG and J. S. SHIPMAN 
(Laboratory jor Electronics Inec.. Boston, Massachusetts, U.S.A.) 


and H. KAUFMAN 
(McGill University, Montreal, P.Q., Canada) 


Received 31 August 1954] 


SUMMARY 


The first few Bennett functions of the second kind, i.e. coefficients A‘) (k) in the 
louble Fourier series expansion of the output from a half-wave square-law rectifier 
a two-frequency input, are tabulated as functions of the input amplitude 


iok. Bennett functions A‘)),(k) of the vth kind for the two-frequency modulation 


" 1 + . > 7. 
product problem for the half-wave vth-law rectifier are defined and recurrence formulae 


mg the A‘ (k) of given kind and different orders, or among those of different 
sinds and given orders, are established for use with the tables. Some brief com- 
ents are made concerning other forms of the problem and of Bennett functions. 


1. Introduction 


THE output from a half-wave square-law rectifier with a two-frequency 


nput of amplitudes P and Pk, 0 < k < 1, can be expressed as a double 
Fourier series in the form 

y(t) > P?A yy (k) he _ A®in(k)cos(w mnt T $200) (1.1) 

mn=0 

where the asterisk indicates that the sum is extended over both signs if 
mn #~ 0 and over only the upper sign if mn = 0 with the single exception 
that the term occurring under the summation sign for m = n = 0 is not 
to be included in the sum. In (1.1), the angular frequencies w ,,,,, and phase 
ingles ¢_,,,,, are particular linear combinations of those associated with the 
input, while the amplitudes are the functions of k defined as 


) all 


mn 


> 


1‘) n(&) — | (cos u+k cos v)? cos mu cos nv dudv, (1.2) 
Js 

where (m, n = 0, 1, 2,...) and R is the region of the square 0 < u,v <7 
in which cosu-+-keosv > 0. For the details we refer to Bennett (1). In 
this paper we give numerical tables of the first sixteen non-zero functions 
Amik) = A®),,(k) as defined by (1.2) for k = 0-02(0-02)1-0. The tables 
given are thus a companion set to the tables given by Sternberg, Shipman, 
ind Thurston (2) for the corresponding problem for the half-wave linear 
rectifier. Each entry in the tables is believed to be accurate to somewhat 
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better than one unit in the last place given and probably most are 
rounded correctly. 


Following the terminology of (2), we refer to the functions A) (k) 
defined by (1.2) as Bennett functions of the second kind. Similarly, if 
the same two-frequency input is applied to a half-wave vth-law rectifier, 
then in place of (1.1) we have for the output the double Fourier series 
expansion 


y(t) = PVAQ(k)+ PY >* AM an(k)cos(w...,f+¢ 3 


mn mn mi ), ( I ad 
mn=0 


and in place of (1.2) we have the Bennett functions of the vth kind defined 


as » 2 
A), (k) — | | (cos u-+-k cos v)” cos mu cos nv dudr, (1.4 
7 J 
where (v = 1, 2,...;m,n = 0, 1, 2,...) and again R is the region of the square 
0< u,v <7 in which cosu+keosv > 0. Recurrence formulae for the 
functions A‘) (k) = A\,,,(k) by means of which the higher-order functions 


of the vth kind may be expressed in terms of the lower-order functions of 
the vth kind are derived. together with a set of recurrence identities between 
the functions A’ )(k) and A\),.(k) of orders m-+-n and m’-+-n’ = m+-n+! 
which make it possible to express similarly the (v+1)th-kind functions in 
terms of those of the vth kind. The tables of the second-kind functions 
given are thus readily extensible both in order and in kind. In the same 
manner, the functions A‘?)(k) tabulated in the present paper may be 
obtained from the functions 4A,,,,(k) = A‘) (k) tabulated by Sternberg, 
Shipman, and Thurston (2) by use of these two sets of recurrence formulae. 
But, owing to round-off errors and loss of order in this process, it is not 
entirely satisfactory; and taking into consideration the great importance 
of the special case of the square-law rectifier, it was decided to tabulate the 
functions A‘?) (‘) separately. 

The construction of the tables is described briefly in section 2 and some 
remarks on the use of the tables, together with the aforementioned 
recurrence identities, are given in section 3. A proof of the latter relations 
is presented in section 4, and in section 5 we conclude with some brie! 
remarks concerning various other forms in which the present problem 01 
Bennett functions may appear, in addition to those treated. 


2. Construction of the tables 

The construction of the tables was carried out on an I.B.M. Card Pr 
grammed Calculator by computational methods similar to those described 
in (2) for the functions A, (k) A‘) (k) and it will not be discussed 1! 


mn mn 


detail. We note, however, for reference the series and formulae employed 
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Thus, firstly, for all values of k the non-zero even-order functions A{;)(‘) 


mn 
yere CO! puted from the rational formulae 
2)/J 2 (2)/}. l y 
LAA 1 +1, AM (hk) i, (2.1) 
Lk. (2.2) 


o 1k) — A2(k) 


For 0-02 0-40 the odd-order functions A‘?)(k) were evaluated from 


13 AX, 7)| (4/3 2 (1/16)k#—...], (2.3) 
A? k/zr)[2+-(1/4)k?+-(1/96)k*+... |, (2.4) 


Ae r)| (4,15 | /3)k?+(3/16)k4—...], (2.5) 
Ae k/z){ (2/3 | /4)k?+-(1/32)k4-+L...], (2.6) 


i‘< k /ar)| (1/2) —(1/24)k?— (1/256)k4—...], (2.7) 
= {\2 ){ (1/12)-+-(1/192)k2+- (3/2560)k4+-...], (2.8) 
Pees (2 r)| —(4/105)+-(1/5)k®?—(5/16)k4*+-...], (2.9) 
ns of AS k/7 (2/15 | /4)k? —(5/32)k*-+...], (2.10) 
wer ASS) (/ k?/ar)| —(1/6 1 /8)k?— (5/256)k4—... |, (2.11) 
1\2 5 /ar)| —(1/12)+-(1/64)k? + (1/512)k4+-... |, (2.12) 
Az (4 /r)[ —(1/96 | /640)k?—(1/2048)k*—...], (2.13) 


t10 
: ic k>/z7)| (1/960 | 2560)k?+- (5/28672)k4+....], (2.14) 
iy | while for 0-42 < kf 0-98 the computation was done from the formulae 


ae, At?) (] 8 /97?)[ (7 +k?) E(k) —(4—4k") K(k)], (2.15) 
va {2)(] 8/972k)[ (1+ 7k?) E(k) — (1+ 2k? — 3k") K(k], (2.16) 
2257?)[ (23 —23k?+ 8k*) E(k)+(—8-+ 12k?—4k*)K(k)], (2.17) 
tet 12 (A 8 /4572k)| (34 7k? — 2k) E(k) + (—34 2k?+k*) K(k), (2.18) 
1i3)(A 8 /457?k?) 24+. 7k®+ 3k4) E(k)+(2—8k?+6h4)K(k)], (2.19) 
412) 8 /225?k3)[ (8 —23k?+-23k4) H(k)- 

+ (—8+ 27k?—24k4+ 15k*)K(k)], (2.20) 
AS?) (] $/110257? 281+ 985k?— 1080k4+- 384k) E(k) 
el (176— 548k? + 564k4— 192k*)K(k)], (2.21) 
{\° 8/157572k)[ (15— 103k2-+- 128k4— 48k*) E(k) 

(— 15+-58k?— 67k4+-24k*)K(k)], (2.22) 
1p AL 8 31572k? 6+ 9k2— 19k4+- 8k*) E(k) 
+-(6—12k?+-10k*—4k*)K(k)], (2.23) 
sed AU 8/31577k)| (8 — 19k?+ 9k4— 6k) E(k) - 
| | (—8-+-23k?— 18k4+ 3k*)K(k)], (2.24) 
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A\®(k) = (8/157572k*)[ (— 48+ 128k?— 103k44+-15k*) E(k) 
+ (48—- 152k?+- 164k4— 60k) K(k)], (2.25) 
Ai (k) = (8/110252k5)| (384— 1080k?+-985k4— 281k*) E(k) 4 
+-(—384-+ 1272k?— 1501k4 +718k8— L05k8)K(k)], (2.26) 
where K(k) and E(k) are the complete elliptic integrals of modulus fk as 


tabulated by Fletcher (3). For k 1-0 the odd-order functions 4°?) (}) 


mn\” 


were computed from the formula 


4‘2) (1-0) es 2.27 
ao [m2—n?]| (m-+-n)?—4]] (m—n)?—4 |x?’ 


where (m+n = 1, 3, 5,...). These series and formulae may be derived by 
methods such as those described in (1) and (2). 

The computation was checked by a number of methods, including 
differencing, and both the basic calculations and the magnitudes of the 
differences indicate that each entry in the tables should be accurate to 
better than 1 x 10-8 units. 


3. Applications of the tables 

For direct applications of the tables linear interpolation should be 
generally good to at least four decimal places everywhere. If greater 
accuracy is desired, a parabolic interpolation formula with second 
differences may be used to obtain accuracy to about six or more decimal 
places, or a five-point Lagrange formula is recommended and should give 
results accurate to about 3 « 10-8 units. 

For applications of the tables to problems requiring values of the higher 
order, or higher-kind functions A‘y)(k‘), we note the recurrence formulae 


(a) (m—n+v+2)AY i n-a (m-+-n 2)AM a n-1—2mkA, 

(b) (m+n-+v)A? —(m—n—v—2)A ,, bie ee i 
(c) (nt+-m-+yv)A!), (n—m 2)A‘Y) 9 —2(n—1)(1/k). 1, —* : 
(d) (n—m+v+2)A™ 1 nad (n+-m—v—2)A™ |, 1—2n(1 Baer, 
for functions of the same kind, and the recurrence formulae 

(a) Ato’? = AY+kA, 

(b) 2mAPEY = (v+1)A™ , ,—(Vt1 AM 12s (3.2 
(c) 2nAMTD = k(v+1)AM, ,—k(v+1)A®, «1 

for functions of different kinds, where in (3.1) we have (v, m, n Te Bee 
m ~ 1in(b):n + 1in(c)) and in (3.2) we have (vy L200, 0 i Re re 
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Qin (b);n 4 Oin (c)). We note in passing that for vy = 1 the relations 
3.1) reduce to the corresponding recurrence relations of (2). 

The application of the relations (3.1) and (3.2) is elementary, and with 
the aid of the tables one may readily compute values of both the higher- 
order and higher-kind functions A‘) (k). It should be noted in this con- 
nexion, however, that owing to round-off errors there will generally be a 
loss of accuracy when k is small in the use of relations (3.1) (c) and (d), while 
for self-evident reasons there will always be a loss of order in the use of 


> 


relations (3.2). 


4, Proof of the recurrence formulae 

To establish the recurrence formulae (3.1) and (3.2) we begin by noting 
that the integrals (1.4) defining the functions A‘) (k) may always be 
written in either of the forms 


» P 


Ainn(k) - R'”(v)cos nv dv, 
ry 
(4.1) 
+) id 9 
A‘ (k) =< | S)’(u)eos mu du- my Lm 
whe re R VY), DO. (et), and [\) are the integrals 
RY (v) | (cosu+kcosv)’ cos mu du, 
(4.2) 
S)'(U) (cosu+k cos v)’ cos nv dv, 
ii 
pot | | (cosu+kcosv)’ cos nv dv|cos mu du, (4.3) 
é ts 
ind where the limits of integration are « = cos-'k, B = cos—1(—kh), 
cos~} keosv), and y(u) cos~!(—(1/k)cos u), 
while in each of these formulae (v = 1, 2,...; m,n = 0, 1, 2,...). Next we note 


by direct substitution and simplification using (4.2) the preliminary identities 


v+1) RY (2 2mk R(v)cos v-+-(m—v—1) RY _,(v) 0, 


m 


(4.4) 
n+v+1)S™ ,(u)+2n(1/k)S'(uw)cos ut+(n—v—1)S)? ,(u) 0, 
where (v, m, 1, 2,...). Multiplying through the first of the relations 


+.4) by unity or cos nv and the second of relations (4.4) by unity or cos mu 
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and integrating the resulting identities, we obtain, after simplification | tabu 


with the aid of (4.1) and (4.3), the further identities (33) 
(a) (m+v+1)A”).,,+mkA\”, 1-4 prob 
OrTe 

| mkA’?, at (m—v 1)A‘ ie 0, * 


(4.5 Ben 
(d) (n+v+1)A?,.,-+n(1/k)AY?. | 


m 


nl exan 


. (v) : (v) P 
+-m(1/k)A' 4. ,+(n—v—1)A®,_, = 0, (4) a 
: 5 . 9 ; . . . rm attre 
where (v a m,n eo 1.2.5 m =~ 0 in (a): n ~ O in (d)). The 
: : : a ; Bent 
recurrence formulae (3.1) (a) for » 1 and (3.1) (d) for m | now follow 
at once by taking n QO orm 0 in the corresponding relation (4.5) while +] 


the proof of the remaining recurrence formulae (3.1) (@), (6), (¢), and (d) may apek. 
be completed readily by eliminating in turn A}, (kh), AS”, (hk), AG, (h), 
and A'y) ,.,(k) from the two relations (4.5) and renumbering the subscripts 
appropriately. aan 
The first of the recurrence formulae (3.2) is a ready consequence of the aR 
definition (1.4) of the functions A‘) (k). To complete the verification of 3. A 
the remaining two recurrence formulae (3.2) (6) and (c), use (4.2) and (44 4. J 
to derive the preliminary identities 5. D 


2m R' 1)(9:) (y-+-1)R% 


’ 100) — 4 1) Ro 1(v), 


(4.6 
2nS*Y(u) = kv+1)S™ (u)—k(v+1)S8™ (a), 
where (v, m, 1, 2,...). Multiplying the first of relations (4.6) by unity 
or cos nv and the second of relations (4.6) by unity or cos mu and integrating 
the resulting identities, we obtain with the aid of (4.1) and (4.3) the remain 
ing recurrence formulae (3.2) (b) and (c). This completes the proofs. 


We note in passing, for reference, the additional identity 


1 1 1 ( 1 v) 1]. 4(v) 1]. ° 
A mn tA‘ l.n | SA m—1'i,n'| kA, .} $k Ay? 2—i!> (4.1 
where (v ee ae 0, 1, 2,...). This may be considered as alternative 


to the identities (3.2). It follows at once from the definition (1.4) of the 
functions A‘) (k) by splitting off one factor from the kernel function of 


mn 


the integral AY, !(k) and rearranging the integrand before performing th 


mn 


integration. 


5. Other aspects of the problem 

In closing this paper it seems worth while to call the reader’s attentio! 
to several variants of form in which the present problem or Bennett 
functions may appear. Thus, to cite two examples, we note, on the ont 
hand, that in treating problems of the type here considered several writers 
have arrived at answers in the form of certain Weber-Schafheitlin integrals 


which seem in several instances to be completely expressible in terms 0! 
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tabulated Bennett functions. See, for instance, formulae (27) and 


33) of Bennett (1). On the other hand, in treating related but different 


nroblems of a statistical nature several authors have derived certain desired 


orrelation functions in forms proportional to what are now tabulated 
Bennett functions without, of course, anticipating the fact. See, for 
73) and (74) of Chapter 3 of Lawson and Uhlenbeck 
4) and formulae (7.14), (7.15), (7.16), and (7.17) of Middleton (5). It is an 


ittractive thought that still further applications of these tables and of 


example, formulae | 


Bennett functions may yet be found.t 


P to the Am Mathematica! Society at their New York meeting, 23 April 
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Tables of the Bennett Functions Aj,(k) 








4616 
"4709 
"4504 
"4901 
*5000 





*4426 SI51 


*4699 s 
"4748 





‘4799 3449 
“4852 7970 





*4905 6024 
“4966 7423 
*5027 1968 


*5089 9450 





"5513 2840 
“5591 4039 
*56071 6067 





5753 8599 
5835 1294 
"5924 3791 
‘6012 5714 
-6102 6663 


*6194 6219 


6288 3937 
6383 9348 
“6451 1952 
“6580 1215 
“6680 6565 
"6782 7385 
*6586 3001 
“6901 26690 
‘7097 5544 
"7205 O619 
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TABLES OF TWO FUNCTIONS REQUIRED IN CERTAIN 
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SUMMARY 


Tables are provided of the functions 
(i) R(X, 7) = eT {(T?—X*)}42 | e-"1,{(w? —X2)}} dw, 


, 
(ii) 2) R(UX,w) dw, 
X 
and an account is given of their derivation. The tables may be useful in the study 
of surges in flow along pipe-lines and in the transmission of electric current along 


cables. 


1. Introduction 
THE problem of pressure surges, due to movements of a valve, in a uniform 
pipe-line conveying viscous liquid was investigated by Binnie (1). Two 
cases were considered: (i) when a valve at the exit is shut instantaneously; 
(ii) when the valve is closed in such a way that the exit velocity falls 
uniformly to zero. Non-dimensional expressions were derived for the 
pressure changes at a section distant 2 from the exit in terms of time / 
and a number & characteristic of the pipe-line and the liquid in it. This 
paper provides tables from which the pressure variation can be readilj 
computed. The analysis is similar to that of surges in electric transmission 
lines, and identical when the leakage is zero (Carslaw and Jaeger, 2); the 
tables may be of use in this case also. It is because of this application to 
unrelated fields that publication in the common ground of a mathematical 
journal has been chosen. 
In case (i) the solution involves the function R(tka, ikt), where 
” 

R(X, T) = eT h{(T?—X?)+-2 | eI {(w?— X24} dw. (1.1 
An expansion in terms of a series of Bessel functions is derived in section? 
which may be used to evaluate R(X, 7) from existing tables of Bess 
functions such as British Association Mathematical Tables, Volume X (3 


(Quart. Journ. Mech. and Applied Math., Vol. VIII, Pt. 4 (1955)] 
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This series reduces to a single term at the most important section which 
is at X 0, where the greatest pressure in the pipe-line is developed. 
In case (ii) the solution involves Q(}kax, $kt), where 


, 
Q(X, T) = 2 | R(X, w) dw. (1.2) 
x 
For this case also, a series of Bessel functions is obtained in section 3, and 
again it reduces to a single term at X 0. 


2. Instantaneous closure of the valve. Case (i) 

The equations of the earlier paper, Binnie (1), will be referred to, and 
they will be indicated by bold type. It is shown in (3.13) that we have 
to evaluate 


t 
R( Aka, $kt) et] St k(t2—a?)N+k | eth fhk(r?—2?)t} dr, (2.1) 
which reduces to R(X, 7') on substituting 
fh Lkt, X tka, w kr. (2.2) 


Thus the number of parameters is diminished by one. With the aid of 
Lommel’s expansion given by Watson (4, p. 140), the integral in (2.1) may 
be put in the form 


1 T 
. : . a u l an” Fs 
ss 2 e* I} (w? . 2 dw 2 S ; 7 > X nT (w) dw. (2.3) 
; j font nm! (2w)" 
X nm 
- 
But I (w) di = U, (ww) I (w)}, (2.4) 
7 bie aT) mw" 
, + \ l )" -_ iad f ' In 
therefore S i. ae oe I, s(w)+1,(w)}]. (2.5) 
9 = =n )u - 
. . 1 )et+iXe . - , 
Now > 2 Un (X)+,(X)} = e®, (2.6) 
fod EL - mit 


n 0 
This identity, which seems to be new, is established in the Appendix, 


section 5. Then we have finally 


can 1)" X* e~* f 1 TV\) VY 
. pS n!2"-1(2n 1) 7” Uh 1 r’) I LAT); -2X 
re ™ 
%e-1 | TST, T)+-1,(7T = {Ty ( D)+1(T)}— Aah?) \ 1,(T) L 
rn 7 : 
= {1,(7')+-1;(7')}—... | —2X. (2.7) 


240T? 
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This result was checked by expanding (2.1), with (2.7) substituted therein, 
into series of polynomials and comparing their sum with the like expansion 
of (3.19), which is the corresponding result in the earlier paper. It was 
during this process that (2.6) was discovered. For X 0, (2.7) reduces to 


one term only. 


3. Exit velocity falling uniformly to zero in time /,. Case (ii) 
It appears from (4.5) that we must evaluate 

z 

Q(x, 7) =2| 

X 


R(X, w) dw. (3.1 


The integral of the first term of R can be directly found from (2.7). That 
of the second term can be established from (2.7) by means of two integrals, 
These are the same as (2.4) with n—1 replacing n, and 


ied u u l uw 


| : (1, (w) dw : f. (wv) : 


Yh yen 


{f,-s(w)+L,-1(w)}, (3.2 
) 2 ‘ « " 


2n—3 w" 
which is obtained by integration by parts. Hence 


a e- ' (a) 2 
ii, (we) +L, (w)} dw | = ~ {Tf -2(w) +L, 1(w)} | 
; 


Th we 2] uw 2n—3 


e* | I 
|i o(w) (2 lI, (2) I (2e) |. (3.3 
n—3) i uw 


ae 2( ») 


There is little to choose between these two forms. We also have the identity 


< (—1)"t+1X" [7 _.(X) 2 7 Yer xy X\ 
— n! 2"(2n—1) X ai alt) nb) | 3) 


n 


which is analogous to (2.6) and is proved in the Appendix, section 6. Hence 





1 
| eee —  (— 1 )n+1 Xen r il i 2 - a 
Sd7 ™ i i ‘ 1 if -(7T)+4,_(T)}I 
a n! 2"-1(2n 1)| 7-2 7 2n—3 i 
+ n= 


xX" A ; X 2n a 


(7) 1(7)}] ! 


(X)+I, (X))] 2X(T—X 


xen (20). 260,(7) +I 7)}) 

X4/1,(7) 

= r 
x6 (A T) 

1207\ T 


Again, at X = 0 only one term remains. 


24 7,(7') 1(7)}] 


/ 


2(7,(7') 1,(7)}] | 2XT—2X3, (35 
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4. Description of the tables 

The tables were computed at the Mathematics Division of the National 
Physical Laboratory, mainly under the direction of Mr. G. F. Miller. The 
table of R(X.7') was compiled by means of the series (2.7) and that of 
Q(X. 7) by quadrature with (3.5) supplying check values and also initial 
values for T’ ~ X. 

The tables give values of R(X, 7') and Q(X, 7’), each for X = 0(0-2)5 
ind 7 = 0(0-2)5(1)20. When 7' < 5, and when X <1, 5 Tl <x 20, 
second differences are prov ided both X-wise and T'-wise for interpolation 
by means of Everett's formula 

bh) = (1—@)f(a)+Of(a+h)+ E5(0)8*f(a)+ By(0)82f(a+-h). (4.1) 
By replacing 5° by (6*—0-1846*) differences were modified where necessary 
to allow for the effect of fourth differences. The formula (4.1) thus always 
suffices to give full accuracy. Differences were rounded off after formation 
nd modification so they cannot be exactly reproduced from the values 
of R and Q given in the tables. 

In those parts of the table where differences are not provided, linear 


interpolation will give an error not exceeding 0-003 in R; the error from 


linear interpolation in Q may be larger, up to 0-15 near X ey 20, 
with a tendency to be smaller for the smaller values of Q@. More accurate 
interpolation is best done in this part of the table by first forming the 
missing differences 

\s an example of the use of the tables, we will evaluate R(0-45, 3-17). 


With the aid of (2.6), (2.8), and (2.2) this example can be shown to corre 


spond to time 88-1 sec. at a point halfway along a pipe-line of diameter 6 in. 


ind length 100,000 ft., which is conveying a liquid of kinematic viscosity 
5:63 10-4 ft.2 see. and in which the wave velocity is, 4,000 ft./sec.) First 
we dete ine R(0-4, 3-17) and R(0-6, 3-17) with 
3°17—3:-0 * 5 
H 0-85. =, 0-15. 
()-? 
B’(0-85) 0-0244., E(0-85) 0-0393. 
‘het R(0-4, 3-17) 2-21939, 3° |- 3132, 
R(0-6, 3-17) 1-89846, 3? +3070, 
the values of 6? being obtained by linear interpolation. Throughout the 
working a fifth decimal is retained as a guard figure. For the X-wise 
nterpolatio1 
0-45 0-4 _ 
6 0-25, 1—é 0-75, 
Q)-2 
E’(0-25) 0-0547, By(0-25) 0-0391, 


whence R(0-45, 3°17) = 2-1362, 
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the guard figure now being rejected. The coefficients £5(@) and E£4(6) were 
taken from Interpolation and Allied Tables (5), pp. 800-801. 


APPENDIX 
5. Derivation of (2.6) 


A proof of the formula (2.6) seems desirable. Consider 


- 8 
r ( -_ yn wes n : . 
3 n! Si - Ly" n_s(X) +1, (X)}. (5.1) 
n=0 
We wish to show that V = eX, and we shall establish this identity by proving that 


dV /dX = V. 


Differentiating (5.1) and using the relations 


dl, (X) x n : dl (X) ? n—1 . 
=—J, ,X)——IJ,(X), = = 1 (X)-+ —— I X), 5.2) 
ax nh )— x Ent) dX SAT yeh = 
we abtain 
dV WD (—1)"4X"{/. Qn 
; — 14 I x Ax 
dX 2. n!2"(2n- Ht X ) n—1(-\) nt-0) i 
: = ])er1Lyr- 1 2 
: 2 19 : I, i(X ). (5.3 
nia” 
n 0 
P : “ n X n ] ‘< 
Now ~ ( 1) Soe e. | = x at 0: (5.4 
Ly tan : s!2° 
n 0 
= ] n Xx” 
* ee | ee oO, 2 l. (5.5 
ya n!2" 
0 
(See Watson (4), p. 141, formulae (7) and (9) with z 7X: the latter formula is true 
only when v > 0, integer or not.) From (5.5) with s 1 we see that the series it 
(5.3) is zero and, since V 1 when X 0, it follows that 
V ex (5.6 


as required, 

In view of the error in Watson’s formula (9) quoted above, it seems worth whil 
to give a short proof of (5.4) and (5.5) for integral s. The generating series for thi 
Bessel functions J,(.X) is 


] NX anil n yn 
whilst exp| “5 S (—)) , 
zt hang YZ" t” 
a 0 
= 4 X 
Tiles m 1)” "n ” 
Thus exp{4.X7¢} - S be nay", (5 
= a _ n gn 
0 


and we may equate coefficients of: positive or zero powers of t on both sides to obtall 
(5.4); we may also equate coefficients of negative powers of t to obtain (5.5). All t! 
series converge absolutely for all values of X and ¢tsince /,(X) ~ (4X)°/s! for large § 
whilst J_,(X) I.(X) when s is an integer. 
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6. Derivation of (3.4) 
This identity may be proved by the method of section 5 applied twice. We write 
W for the left side of (3.4) and differentiate, obtaining 


dV = 1n41X"~ 2 
Ss T,_o(X)+Iy_.(X)}4 


dX Ly nt 2"(2n—1)L2n—3*" =“ 


) 


2n Pe 4 2 R 3 
oeone Steet _ 25) w—2, (6.1) 


as 5 Ms 








‘ ] w+1lyn 1 . - 
shew See > ke ih Ee. (6.2) 
Ly n!2"(2n—3) . 

dZ =m (—])nt+iXn-2 : ; uf 
Again z+ con i AD an O (6.3) 

dX y- n!2” 
the series being zero by (5.5) with s 2: therefore 

Z Xe x, 

and W (w—AX Je x, (6.4) 


n which A and p are constants. But at X 0, W 1 and dW/dX 4/3, whence 
3.4) follows 
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Tabli of Q(X, r) 
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R 83, R 8 | 8 R 8 | 8h 
X 16 18 2-0 

o'2019 | —15 276 
*2739 14 271 oO 1653 6 245 

| 

0°3445 | 13 267 | 0°2278 6 |+ 242 | 0°1353 o |+217 
4138 | 12 263 2897 6 240 *1895 I 216 
*4520 II 259 3509 6 237 2435 I 214 
*5490 II 255 4114 6 234 2974 2 213 
“6149 10 251 4714 6 232 “3511 2 212 
0°6797 10 248 | 0°5307 6 |+ 229 | 074045 3 |+210 
*7430 9 244 5893 6 227 °4577 3 209 
+8066 9 241 *6473 6 225 "5105 3 207 
*8686 9 238 "7047 6 222 +5631 3 206 
*9295 8 235 “76015 6 | 220 “60153 4 204 
o'QgoI 232 | 0°8177 6 |+218 | 06671 4 |-+ 203 
1°0496 8 8733 6 215 7185 4 201 
1'1084 7 "9283 6 213 7696 rt 200 
1°1664 7 0°9525 6 211 8203 4 198 
1°2237 7 1°0307 5 209 8706 4 197 
1°2502 7 218 1*Ogoo 5 207 | 0°9205 4 |+195 
R 383, R 8} by R 8} + R 

za 2°4 26 28 30 

yo ee 5 Lainie = a 
907 41-8 +. 768 * 

1310 6 169 743 | +10 | -+- 147 * 

1719 5 170 1059 9 149 605 12 t- 128 

+ 

2134 4 |) 171 11444 7 | 15! go 10 |-+- 131 498 | +12 

2552 3 171 | 1807 6 152 | 1213 9 133 753 II 

2974 2 172 | 2175 5 153 | 1529 7 135 |} 1018 9 
3398 2 172 | 2548 4 154 | 1552 6 137 | 1294 8 
3824 I 172 | 2926 3 155 | 2152 5 138 | 1577 7 
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di 
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Oy 

112 * 
115 | 408 
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124 | 1342 
125 1595 
127 18600 
125 2129 
129 | 2403 








* In units of the 4th decimal. 
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